We claim: 


1 1. A method of automatically analyzing at least one seedling germinated from at least one 

2 seed, comprising the steps of:\ 

3 (a) capturing a digital image of the at least one seedling; 

4 (b) identifying the at leasKone seedling in the captured digital image; 

5 (c) determining a primary path of the at least one seedling; 

6 (d) determining at least one Value from the primary path of the at least one seedling; 

7 and \ 

8 (e) determining a seed vigor inaex .from at least the at least one value determined 

9 from the primary path of the at least one seedling. 

1 2. The method of automatically analyzing at least one seedling according to claim 1 : 

2 (a) wherein said step of determining aNeast one value from the primary path of the at 

3 least one seedling comprises the step of determining a value corresponding to an 

4 overall length of the at least one seedling\from the primary path of the at least one 

5 seedling; and \ 

6 (b) wherein said step of determining a seed vigor index from at least the at least one 

7 value determined from the primary path of the arJeast one seedling comprises the step 

8 of determining a seed vigor index from at least thk value corresponding to the overall 

9 length of the at least one seedling. \ 

1 3. The method of automatically analyzing at least one seedlmg according to claim 1 further 

2 comprising the step of determining a separation point between\he hypocotyl of the at least 

3 one seedling and the radicle of the at least one seedling; and: \ 

4 (a) wherein said step of determining at least one value from \he primary path of the at 

5 least one seedling comprises the step of determining a value corresponding to the 

6 length of at least one of the hypocotyl of the at least one seedlin&and the radicle of the 

7 at least one seedling; and \ 
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8 (b) wherein said^tep of determining a seed vigor index from at least the at least one 

9 value determined from the primary path of the at least one seedling comprises the step 

10 of determining a seecl vigor index from at least the value corresponding to the length 

11 of at least one of the hypocotyl of the at least one seedling and the radicle of the at 

12 least one seedling. \^ 

1 4. The method of automatically analyzing at least one seedling according to claim 1 further 

2 comprising the step of determining a separation point between the hypocotyl of the at least 

3 one seedling and the radicle of the at least one seedling; and: 

4 (a) wherein said step of determining at least one value from the primary path of the at 

5 least one seedling comprises the step of determining a hypocotyl length value 

6 corresponding to the length of the hypo&otyl of the at least one seedling and a radicle 

7 length value corresponding to the length o^ the radicle of the at least one seedling; and 

8 (b) wherein said step of determining a seedNvigor index from at least the at least one 

9 value determined from the primary path of the 'at least one seedling comprises the step 

10 of determining a seed vigor index from at least\the hypocotyl length value and the 

1 1 radicle length value. 

1 5. The method of automatically analyzing at least one seedling according to claim 1 

2 wherein said step of determining a primary path of the at least o^e seedling comprises the step 

3 of determining a locus of pixels, the locus of pixels correspondingsto the primary path of the 

4 at least one seedling and the locus of pixels being narrower in width\than the width of the at 

5 least one seedling in the digital image of the at least one seedling. \^ 

1 6. The method of automatically analyzing at least one seedling according to claim 1 

2 wherein said step of determining a primary path of the at least one seedling comprises the step 

3 of determining a locus of pixels, the locus of pixels corresponding to the primary path of the 

4 at least one seedling and the locus of pixels being a predetermined number of pixels in width. 
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1 7. The method of automatically analyzing at least one seedling according to claim 1 

2 wherein said step of determining a primary path of the at least one seedling comprises the step 

3 of determining a locus of pixel^the locus of pixels corresponding to the primary path of the 

4 at least one seedling and the locusW pixels being one pixel in width. 

1 8. The method of automatically analyzing at least one seedling according to claim 1 further 

2 comprising the step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling, and \ 

4 (a) wherein said step of determining a primary path of the at least one seedling 

5 comprises the step of determining a primary path for each of the separately identified 

6 overlapped seedlings; \ 

7 (b) wherein said step of determining at least one value from the primary path of the at 

8 least one seedling comprises the step of determining from the primary path for each of 

9 the separately identified overlapped seedlinks a value corresponding to an overall 

10 length of that separately identified overlapped sfeedling; and 

1 1 (c) wherein said step of determining a seed vigo\ index from at least the at least one 

12 value determined from the primary path of the at least one seedling comprises the step 

13 of determining a seed vigor index from at least the plurality of values determined in 

14 step (b). \ 

1 9. The method of automatically analyzing at least one seedling according to claim 1 further 

2 comprising the step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling, and \ 

4 (a) wherein said step of determining a primary path of the at least one seedling 

5 comprises the step of determining a primary path for each of the Separately identified 

6 overlapped seedlings; \ 
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7 (b) wherein said step of determining at least one value from the primary path of the at 

8 least one seedling comprises the step of determining from the primary path for each of 

9 the separately identified overlapped seedlings a value corresponding to the length of at 

10 least one of the hypocotyl of that separately identified overlapped seedling and the 

1 1 radicle of that separately identified overlapped seedling; and 

12 (c) wherein said step of determimng a seed vigor index from at least the at least one 

13 value determined from the primary\path of the at least one seedling comprises the step 

14 of determining a seed vigor index mom at least the plurality of values determined in 

15 step (b). \ 

1 10. The method of automatically analyzing at least one seedling according to claim 1 further 

2 comprising the step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling, and \ 

4 (a) wherein said step of determining a primary path of the at least one seedling 

5 comprises the step of determining a primary path for each of the separately identified 

6 overlapped seedlings; \ 

7 (b) wherein said step of determining at least one value Nfrom the primary path of the at 

8 least one seedling comprises the step of determining froni the primary path for each of 

9 the separately identified overlapped seedlings a hypocotyl length value corresponding 

10 to the length of the hypocotyl of that separately identified oVerlapped seedling and a 

1 1 radicle length value corresponding to the length of the radkle of that separately 

12 identified overlapped seedling; and \ 

13 (c) wherein said step of determining a seed vigor index from at lelast the at least one 

14 value determined from the primary path of the at least one seedling comprises the step 

15 of determining a seed vigor index from at least the plurality of hybocotyl length 

16 values determined in step (b) and at least the plurality of radicle length values 

17 determined in step (b). \ 
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1 11. The method of automatically analyzing at least one seedling according to claim 8 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedliiig comprises evaluating an energy function. 

1 12. The method of automaticd^ analyzing at least one seedling according to claim 9 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function. 

1 13. The method of automatically analyzing at least one seedling according to claim 10 

2 wherein said step of separately identifying aNplurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function. 

1 14. The method of automatically analyzing at least one seedling according to claim 8 

2 wherein said step of separately identifying a plurality^of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings. 

1 15. The method of automatically analyzing at least one seedling according to claim 9 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 


r fW 

4 configurations of at least partial primary paths of overlapped seedlingk 


3 image of the at least one seedling comprises evaluating an energy function based on proposed 


1 16. The method of automatically analyzing at least one seedling according to claim 10 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings. 

1 17. The method of automatically analyzing at least one seedling according ^to^ claim 8 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the^digital 
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3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings using the following 

5 heuristics: primary paths do not make unnaturally sharp turns and seedling edges should 

6 be used as much as possible. 

1 18. The method of automatically Analyzing at least one seedling according to claim 9 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings using the following 

5 heuristics: primary paths do not make unnaturally sharp turns and seedling edges should 

6 be used as much as possible. 

1 19. The method of automatically analyzing ak least one seedling according to claim 10 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings using the following 

5 heuristics: primary paths do not make unnaturally sharp turns and seedling edges should 

6 be used as much as possible. 

1 20. The method of automatically analyzing at least one\seedling according to claim 8 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energ^function based on proposed 

4 configurations of at least partial primary paths of overlapped seedlings using the following 

5 heuristics: primary paths should not make unnaturally sharp turns, seedling edges should 

6 be used as much as possible, and all primary axes should hav v e a hypocotyl/radicle 

7 separation point. 
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1 21. The method of automatically analyzing at least one seedling according to claim 9 

2 wherein said step of separatelyNidentifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedlingWnprises evaluating an energy function based on proposed 

4 configurations of at least partial pimiary paths of overlapped seedlings using the following 

5 heuristics: primary paths should not\make unnaturally sharp turns, seedling edges should 

6 be used as much as possible, and all primary axes should have a hypocotyl/radicle 

7 separation point. 

1 22. The method of automatically analyzing at least one seedling according to claim 10 

2 wherein said step of separately identifying a plurality of overlapped seedlings in the digital 

3 image of the at least one seedling comprises evaluating an energy function based on proposed 

4 configurations of at least partial primary paths of Overlapped seedlings using the following 

5 heuristics: primary paths should not make unnaturally sharp turns, seedling edges should 

6 be used as much as possible, and all primary axes should have a hypocotyl/radicle 

7 separation point. 


1 23. The method of automatically analyzing at least one seedling according to claim 1 further 

2 comprising the steps of determining a first locus of points indicating the hypocotyl of the at 

3 least one seedling, determining a second locus of points indicating the radicle of the at least 

4 one seedling, overlaying the first and second loci over an image ofNjie seedlings to generate a 

5 composite image, and displaying the composite image. 

1 24. The method of automatically analyzing at least one seedling according to claim 23 

2 wherein said step of displaying the composite image comprises the step^ of displaying the 

3 composite image on a video display terminal. 
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1 25. The method ofWtomatically analyzing at least one seedling according to claim 23 

2 wherein said step of displaying the composite image comprises the step of printing the image 

3 on a printer or plotter. \ 

1'" 26. A method of automatically analyzing at least one seedling germinated from at least one 

2 seed, comprising the steps of: 

3 (a) capturing a digital image of the at least one seedling; 

4 (b) determining a first locus of points indicating the hypocotyl of the at least one 

5 seedling; 

6 (c) determining a second locus of points indicating the radicle of the at least one 

7 seedling; 

8 (d) overlaying the first and second loci over an image of the seedlings to generate a 

9 composite image; and 

10 (e) displaying the composite image. 

A^27. A method of analyzing at least one seedling germinated from at least one seed, 

2 comprising the steps of: 

3 (a) placing a growing medium in a shallow container; 

4 (b) wetting the growing medium; 

5 ^ (c) placing the at least one seed onto the growing medium; 

6 (d) germinating the at least one seed with the shallow container at an angle with 

7 respect to the vertical that is less than about 10°; 

8 (e) capturing a digital image of the at least one seedling; and 

9 (f) analyzing the captured digital image of the germinated seedling. 


\( 28^ 


The method of analyzing at least one seedling according to claim 1 wherein said step of 
2 germinating the at least one seed with the shallow container at an angle with respect to the 
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3 vertical that is less than about 10° comprises the step of positioning the shallow container 

4 vertically. 

1 29. The method of analyzing at least one seedling according to claim 1 wherein said step of 

2 capturing a digital image of the at least one seedling comprises capturing an image of the at 

3 least one seedling using a scanner having a scanner surface and positioned with its scanner 

4 surface oriented at least 90° from the horizontal. 

1 30. The method of analyzing at least one seedling according to claim 1 wherein said step of 

2 capturing a digital image of the at least one seedling comprises capturing an image of the at 

3 least one seedling using a scanner having a scanner surface and positioned with its scanner 

4 surface substantially inverted so that it captures an image of at least one seedling positioned 

5 beneath the scanner. 
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// SeedlingAnalyzer . cpp : implementation of the SeedlingAnalyzer class. 
// 

////////////////////////////////////////////////////////////////////// 

#include "stdafx.h" 

#inc lude " Progeny . h " 

#include "SeedlingAnalyzer . h n 

# include "ProgenyDoc .h" 

#include "ColorConversion. h" 

#include "const. h" 

#include <math.h> 

#include <f stream. h> 

#ifdef _DEBUG 
#undef THIS_FILE 

static char THIS_FILE [] = FILE ; 

#define new DEBUG_NEW 
#endif 

4fdefine GETRAND (rand ()/ (double) (RAND MAX+1) ) 

;,y — 

$//////////////////////////////////////////////////////////// ///////// 

Construction/Destruction 
77//////////////////////////////////////////////////////////////////// 

^SeedlingAnalyzer : : -SeedlingAnalyzer ( ) 

rip 

Q 

r^/ Returns information needed to measure radicle and hypocotyl lengths on 
s.^/ each seedling detected from a color image of seedlings . 

rgtd : : vector<SeedlingInf o *> SeedlingAnalyzer :: ProcessSeedlinglmage (ImagelO* image) 

>:i 

// Minimum pixel area for a blob to be qualified as a seedling blob 

const int c_minblobsize = 100; 

// kernel size used for median filtering 

const int c_medkernelsize = 3; 

// Minimum skeleton length to be qualified as a seedling blob 
const int c_minseglen = 15; 

// Minimum size for a blob to be considered as a cotyledon 
const int c_mincotsize = 10; 

// Minimum size for a blob to be considered as a seed coat 
const int c_mincoatsize =10; 

// Seedling separation parameters 

const float c_angleW = 200. Of; 

const float c_initialTemperature = 2000. Of; 

const float c_lengthW = O.Of; 

const int c_maxIteration = 50000; 

const float c_temperatureConst = l.Of; 

const float c_unusedW = 10. Of; 

const float c_minEdgeLength = 10. Of; 

const float c_separationW = 2 00. Of; 

// Data structure to hold extract measurements of seedlings 
std: : vector<Seedling!nfo *> seedlinginf o; 
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// origimg is the original RGB image of seedlings. 
Exlmage* origimg = (Exlmage* ) image ; 

// Create binary image of cotyledons (leaves) 

Exlmage* cotimg = origimg- >Threshold (2 00 , 255, 200, 255, 0, 200); 
// Create binary image of seed coats 

Exlmage* coatimg = origimg- >Threshold (60 , 255, 0, 255, 0, 100); 
// Binary threshold based on red channel 

//To threshold, find the peak intensity on the red channel by FindPeakValue (0) . 

// The peak intensity + 40 is a good lower bound for thresholding. 

Exlmage* binimg = origimg- >Threshold (origimg- >FindPeakValue (0) +40 , 255, 0, 255, 0, 


// Remove blobs of size less than c_minblobsize 

BinaryObjects* binobjs = binimg- >GetObj ects () ; 

delete binimg; 
Q binobjs->FilterOut (c_minblobsize) ; 
•:.g Exlmage* cleanbinimg = binobj s->CreateImage ( ) ; 
I.n delete binobjs; 

M> // Perform median filtering to remove noise and smooth out edges 

l^j Exlmage* smoothimg = cleanbinimg- >MedianFilter (c_medkernelsize) ; 
\U BinaryObjects* smoothobjs = smoothimg- >GetObj ects () ; 

Q std: : vector<0bjectlnf o*> blobobjinfo = smoothobjs->GetObjectInfo () ; 

en 

i»3 // 4. Perform thinning (skeletonization) 

■• n 

^ Exlmage* thinimg = smoothimg- >Thin () ; 

^ delete smoothimg; 

// Throw out skeletons of length less than c_minsegln 
BinaryObjects* binskel = thinimg- >GetObj ects () ; 
delete thinimg; 

binskel->FilterOut (c_minseglen) ; 

// objinfo holds information for all detected objects in the 
// skeletonized image. 

std: :vector<0bjectlnfo*> objinfo = binskel->GetObjectInf o ( ) ; 

// Since the extracted skeleton contains root hairs and other 

// noise, use the shortest path algorithm to find the primary 

// axis. Also, the second junction in the path is the hypocotyl 

// separation, if it is not the end junction (i.e., not the bottom) . 

ShortestPathFinder spf; 

ConnectivityGraph<Edge *>* jg; // Junction graph. 

bool added = false; // True if seedlinginfo was updated. 

// Perform the following loop for each seedling blob 
for(i=0; i<obj info . size () ; i++) 

{ 
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// Compute the junction graph for the seedling blob 
std: :vector<Junction*> jlist; 

jg = binskel->CreateJunctionGraph (obj inf o [i] , jlist); 

// A seedling blob has to have more than 2 junctions 
// to be considered for further processing 
if (jlist. sizeO > 1) 

{ 

// Find the cotyledon junction. 

// Find cotyledons and leaves in the current blob. 

Exlmage* cotblobimg = cotimg->Crop (max (0 , obj inf o [i] ->xmin-8) , 
max (0, obj info [i] ->ymin-8) , 

min (origimg->GetWidth () -1, objinfo [i] ->xmax+8) , min (origimg- 
>GetHeight () -1, obj info [i] ->ymax+8) ) ; 


Exlmage* coatblobimg = coatimg->Crop (max (0 , obj inf o [i] ->xmin-8) , 
.max (0, obj info [i] ->ymin-8) , 

^ min (origimg->GetWidth () -1, objinfo [i] ->xmax+8) , min(origimg- 

: =js>GetHeight () -1, obj info [i] ->ymax+8) ) ; 

m 

// Throw out blobs of size less than c__mincotsize 

~ // cotinfo contains information for all extracted cotyledon 

; 5 h // blobs greater than or equal to c_mincotsize 

™ BinaryObjects* cotobjs = cotblobimg- >GetObj ects () ; 

delete cotblobimg; 
U cotobj s->FilterOut (c_mincotsize) ; 

l!fi std: : vector<0bjectlnf o*> cotinfo = cotobj s->GetObject Inf o () ; 

v3 // Throw out blobs of size less than cjmincoatsize 

Q // coatinfo contains information for all extracted seed coat 

□ // blobs greater than or equal to c_mincoatsize 

BinaryObjects* coatobjs = coatblobimg- >GetObj ects () ; 

delete coatblobimg; 

coatobjs ->FilterOut (c_mincoatsize) ; 

std: :vector<0bjectlnfo*> coatinfo = coatobj s->GetObjectInfo ( ) ; 

// coatindx holds id of seed coat junctions 
std: :vector<int> coatindx; 


// Find the junction (in junction graph) that corresponds to seed coats 
for(int c=0; c<coatinf o . size ( ) ; C++) 

{ 

// Since objects were extracted from a cropped image, add offset 
int coatx = max (0 , obj info [i] ->xmin-8) + (coatinfo [c] ->xmax + 

coatinfo [c] ->xmin) /2 ; 

int coaty = max(0,objinfo [i] ->ymin-8) + (coatinfo [c] ->ymax + 

coatinfo [c] ->ymin) /2 ; 


float mindist = 1000000. Of; 
int minindx = -1; 

// Loop through each junction in the junction graph to 
// find the closest junction to the current seed coat 
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(coaty - junc->m_y) * (coaty 
from the seed coat blob 


for(int j=0; j <j list . size () ; j++) 

{ 

Junction* junc - jlist[j]; 

int dist = (coatx - junc->m_ x) * (coatx - junc->m_x) + 
junc->m_y) ; 

// the closest junction must be less than 10 pixels away 


Q 

in 


ru 


if ((dist < mindist) && (dist < 10*10)) 
{ 

mindist = dist; 
minindx = j ; 

} 

} 

// if such a seed coat is found, add it to coatindx 
if (minindx >= 0) 

{ 

coat indx.push_back (minindx) ; 

} 


// cotindx holds id of seed coat junctions, 
std: :vector<int> cotindx; 

// Find the junction that corresponds to cotyledons. 
for(c=0; c<cotinfo. size () ; C++) 


{ 


E iBotinf o [c] - >xmin) 1 2 ; 
Cgotinfo[c] ->ymin)/2; 


cotindx. push_back (-1) ; 

// Since objects were extracted from a cropped image, add offset, 
int cotx = max (0 , obj info [i] ->xmin-8) + (cotinfo [c] ->xmax + 

int coty = max (0 , obj info [i] ->ymin-8) + (cotinfo [c] ->ymax + 


TRACE3 ("blob %d: cotx=%d, coty=%d\n", i, cotx, coty); 
float mindist = 1000000. Of; 

// Loop through each junction in the junction graph to 
// find the closest junction to the current cotyledon 
for(int j=0; j <j list . size () ; j++) 

{ 

Junction* junc = jlist[j]; 

int dist = (cotx - junc->m_x) * (cotx - junc->m_x) + (coty 
- junc->m_y) * (coty - junc->m_y) ; 

if (dist < mindist) 

{ 

mindist = dist; 
cotindx [c] = j; 

} 

} 

} 


cotyledon 


// Each seed coat that is close enough to a cotyledon is merged to the 
for(int cot=0; cot<cotindx . size ( ) ; cot++) 
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for(int coat=0; coat<coatindx. size ( ) ; coat++) 

{ 

// if junctions are within 10 pixel distance, then drop the 

coat junction. 

float sqdist = pow (j list [cotindx [cot] ] ->m_x - 
j list [coatindx [coat] ] ->m_x, 2) + 

pow(jlist [cotindx [cot] ] ->m_y - jlist [coatindx [coat] ] - 

>m_y, 2) ; 


£3 
5 


if(sqdist < 10.0f*10.0f) 

{ 

coatindx. erase (coatindx.beginO +coat) ; 
coat-- ; 

} 


} 


// Free memory 

for(int obj=0; obj <cotinf o . size ( ) ; obj++) 
£3 delete cotinf o [obj ] ; 

! =n // Free memory 

i;p for(obj=0; obj<coatinfo. size () ; obj++) 

M delete coatinf o [obj ] ; 

=3 

// primaryAxes holds a primary axis for each detected seedling 
\U II in the current seedling blob 

JB std: :vector<std: :vector<int> > primaryAxes; 

i'fj // If no cotyledon/ seed coat was detected, 

i«3 // assume the seedling blob contains only one seedling 

if ( (cotindx. size () == 0) && (coatindx. size () == 0)) 

{ 

// Find the longest shortest path to find the primary axis 
float maxdist = O.Of; 
int maxindxl = -1, maxindx2 = -1; 
for(int k=0; k< jg- >GetSize ( ) ; k++) 

{ 

std: : vector<f loat> dist = spf . FindShortestPathDistance (* jg, 
for(int m=0; m<dist . size ( ) ; m++) 

{ 

if (maxdist < dist [m] ) 

{ 

maxdist = dist [m] ; 
maxindxl = k; 
maxindx2 - m; 

} 

} 

} 

ASSERT ( (maxindxl != -1) && (maxindx2 != -1) ) ; 

// Since no information is given on which end junction belongs 
// to a cotyledon, find out which index is top/bottom. 
// Assume the top junction belongs to a cotyledon, 
int yl = jlist [maxindxl] ->m_y; 
int y2 = j list [maxindx2 ] - >m_y; 
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if(yl > y2) 
{ 

int temp = maxindxl; 
maxindxl = maxindx2 ; 
maxindx2 = temp; 

} 

primaryAxes .push^back (spf . FindShortestPath (* jg, maxindxl , 

} 

// If a cotyledon was detected but no seed coat, find the shortest 
// path from the cotyledon junction to every other junciton, 
// and take the longest path as the primary axis 
else if ( (cotindx. size () == 1) && (coatindx. size () ==0) ) 

{ 

std: :vector<float> dist = spf . FindShortestPathDistance (*jg, 


int maxindx = -1; 
float maxdist = -l.Of; 

for (int m=0; m<dist . size () ; m++) 

{ 

if (maxdist < dist [m] ) 

{ 

maxdist = dist [m] ; 
maxindx - m; 

} 

} 

ASSERT (maxindx >= 0) ; 

int yl = j list [cotindx [0] ] ->m_y; 
int y2 = j list [maxindx] - >m_y; 

int indxl, indx2 ; 

if (yl > y2) 
{ 


} 

else 

{ 


indxl = maxindx; 
indx2 = cotindx [0] ; 


indxl = cotindx [0] ; 
indx2 = maxindx; 


} 

primaryAxes .push_back (spf . FindShortestPath (*jg, indxl, indx2)) ; 

} 

// If there is one seed coat and no cotyledon, assume one seedling 
// within the blob 

else if ( (coatindx. size () == 1) && (cotindx . size () ==0) ) 

{ 

// Find the longest shortest path. 

std: :vector<float> dist = spf . FindShortestPathDistance (*jg, 
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int maxindx = -1; 
float maxdist = -l.Of; 

for(int m=0; m<dist . size ( ) ; m++) 

{ 

if (maxdist < dist [m] ) 

{ 

maxdist = dist [m] ; 
maxindx = m; 

} 

} 

ASSERT (maxindx >= 0) ; 

int yl = j list [coatindx [0] ] ->m_y; 
int y2 - j list [maxindx] ->m_y; 

int indxl, indx2; 

if (yl > y2) 
{ 

indxl = maxindx; 
indx2 = coatindx [0] ; 


■■■-r else 

w { 

,u indxl = coatindx [0] ; 

; = indx2 = maxindx; 

2 > 

La 3 

□ primaryAxes .push back (spf . FindShortestPath (* jg, indxl, indx2) ) ; 

.0 } 

Q // This is the case where there are multiple cotyledons in the blob 

f3 else if (coatindx. size () + cotindx. size () > 1) 

{ 

std: :vector<int> start juncs ; 
for(int c=0; c<coatindx . size ( ) ; C++) 

start juncs ,push_back (coatindx [c] ) ; 
for(c=0; c<cotindx. size {) ; C++) 

start juncs .push_back (cotindx [c] ) ; 
// Obtain primary axis for each seedling in the blob 
primaryAxes = SeparateSeedlings ( jg, startjuncs, c_maxIteration, 
c_lengthw, c_angleW, c_unusedW, c^initialTemperature, c_temperatureConst , 
c_minEdgeLength, c_separationW) ; 

} 

// Perform hypocotyl/radicle separation on each seedling found in the 

blob 

// The first junction encountered while traversing the primary path 

from the 

// cotyledon junction that separates the seedling into hypocotyl and 

radicle such that 

// hypocotyl : radicle length ratio is no less than 0.15 
// is the separation point 

for(int s=0; s<primaryAxes . size ( ) ; s++) 

{ 


22727/04060 


7 


Listing 1 


if (primaryAxes [s] . size () > 1) 
{ 

int k = 1; 

Seedlinglnfo* sinfo = new Seedlinglnf o; 

sinf o->hypocotyl_length = j g->Get Edge (primaryAxes [s] [k-1] , 

primaryAxes [s] [k] ) ->m_length; 

sinfo- >radicle_length = 0; 

for(int j=l; j <primaryAxes [s] . size ( ) -1 ; j++) 

{ 

sinfo- >radicle_length += jg- 
>Get Edge (primaryAxes [s] [j] , primaryAxes [s] [j+1] ) ->m_length; 

} 

// Check to see if hypocotyl/radicle ratio is not so 

extreme . 

// If so, extend hypocotyl and shorten radicle. 
k++; 

; =f while ( ( (float) sinf o->hypocotyl_length / sinfo- 

! :-Pradicle length < 0.15) && (primaryAxes [s] . size ( ) > k) ) 

m ~ { 

V3 sinf o->hypocotyl_length += jg- 

j=£GetEdge (primaryAxes [s] [k-1] , primaryAxes [s] [k] ) ->m_length; 
:I p sinf o->radicle_length -= jg- 

^WGetEdge (primaryAxes [s] [k-1] , primaryAxes [s] [k] ) ->m_length; 

} 


primaryAxes [s] [m+1] ) 


sinf o->hyporad_separat ion = k - 1; 
sinf o->primary_axis = primaryAxes [s] ; 
sinfo- >junction_graph = jg; 

// Compute the bounding box for the seedling 

int xmin = 1000000, xmax = -1, ymin = 1000000, ymax = -1; 
for(int m=0; m<primaryAxes [s] . size ( ) -1 ; m++) 

{ 

Edge* edge = j g->Get Edge (primaryAxes [s] [m] , 

if (edge->m_xmin < xmin) 

xmin = edge->m_xmin; 
i f ( edge - >m_xmax > xmax ) 

xmax = edge - >m_xmax ; 
if (edge->m_ymin < ymin) 

ymin = edge - >m_ymin ; 
if (edge->m_ymax > ymax) 

ymax = edge - >m_ymax ; 

} 

sinf o->toplef t .x = xmin; 
sinfo->toplef t .y = ymin; 
sinfo- >bottomright .x = xmax; 
sinf o->bottomright .y = ymax; 

seedlinginfo .push_back (sinf o) ; 
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added = true; 

} 

} // if there is one cotyledon in the blob 
} // if there are more than one junctions 
if ( 'added) 

delete jg; 

else 

added = false; 

} 

// In case cotyledon/hypocotyl separation was undetected, 
// assume separation point at the mean ratio. 

for(int i=0; i<seedlinginf o . size ( ) ; i++) 

{ 

if (seedlinginf o [i] ->hyporad_separation == 0) 

{ 

seedlinginfo [i] ->useAverageSeparation = true; 

} 

Q } 

in / / Free memory . 

M delete cotimg; 

;a p delete coatimg; 

! s y return seedlinginfo; 

"L 

Separate seedlings using simulated annealing. 
:=atd: :vector<std: :vector<int> > SeedlingAnalyzer : : SeparateSeedlings (Connect ivityGraph<Edge 

jg ra P n / std: :vector<int> start junc, 
pint loopmax, float lengthW, float angleW, float unusedW, float temperature, float 
^emperatureConst , float minEdgeLength, float separationW) 

: t 

// start contains the ID of junctions that are preassigned to seedlings. 
// start. size () is the number of seedlings assigned initially. 

// Keep track of the path for each seedling, 
std: :vector<std: :vector<int> > paths; 

// Keep track of the length for the paths, 
std: :vector<f loat> pathlength; 

// Keep track of the last "valid" end angle taken for each seedling, 
std: :vector<std: :vector<f loat> > endAngle; 

// Keep track of hypocotyl/radicle separation 
std: :vector<int> RHseparation; 
std: :vector<int> RHseparation2 ; 
std: : vector<f loat> separationlen; 

// Initialize data structures by starting out with a single 
// junction for each seedling 
for(int i=0; i<start junc . size () ; i++) 

{ 

std: :vector<int> path; 
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# 


paths .push_back (path) ; 

paths [i] .pushjback (start junc [i] ) ; 

std: : vector<f loat> list; 
endAngle.push_back (list) ; 

pathlength.push_back(0 .Of ) ; 

RHseparation.push_back (-1) ; // Initially, paths don't have RH separation 
RHseparation2 ,push_back (-1) ; 
separationlen .push_back (0 . Of ) ; 


// Keep track of which edge is occupied 

std: :vector<std: : set<int> > edgeOccupation; 
for(i=0; i<j graph- >GetNumEdges ( ) ; i++) 

{ 

std: :set<int> edgeset; 
edgeOccupat ion. push_back (edgeset) ; 

■5 } 


not . 


// Set initial temperature for annealing, 
float energy = O.Of; 


-S-=7 

i s U // Compute the energy of the configuration. 

;: // Since no edges are occupied by seedlings yet, 

□ // the energy is the sum of penalties for 

rf: // unoccupied edges. 

i'3 for(i = 0; i< j graph->GetSize ( ) ; i++) 

{ 

for(int j=i+l; j < jgraph->GetSize ( ) ; j++) 

{ 

if ( jgraph->IsEdge (i, j)) 

energy += unusedW * j graph ->Get Edge (i, j ) ->m_length; 

} 

} 

// add energy for not having hypocotyl/radicle separation, 
energy += separationW * start junc . size () ; 

for(i=0; i<loopmax; i++) 

{ 

// Choose which seedling to update. 

int seedID = (int) (GETRAND* (double) start junc . size ()) ; 
int endjunc = paths [seedID] . back () ; 


bool changeAngle = false; // Whether to change the last "valid" angle or 

// Extend end or remove end. 

// Choose an edge by throwing a die. 

// Copy all neighboring junctions to neighbors. 

std: : vector<int> neighbors; 

for(int j=0; j < jgraph->GetSize ( ) ; j++) 
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{ 

if (j graph- >IsEdge (end j uric, j ) ) 
neighbors ,push_back (j ) ; 

} 

int choice = neighbors [ (int ) (GETRAND* (double) neighbors . size ())] ; 

// Try to remove edge if choice is the junction one before the 
// current junction. 

if ( (paths [seedID] . size () > 1) && (choice == * (paths [seedID] .end() -2) ) ) 

{ 

// Remove edge compute delta energy. 

// * > . . . * >* 

// start choice end 

float deltaEnergy; 

float edgeLength = j graph- >GetEdge (choice, endjunc) - >m_length; 

// Removing hypo/rad separation increases energy 
if (RHseparation2 [seedID] == choice) 

□ { 

=9 deltaEnergy += separationW; 

j\ } 

n 

■■F // If the edge is sufficiently long, subtract energy for the angle 

.-J if (edgeLength > minEdgeLength) 

v { 

if (endAngle [seedID] . size ( ) > 1) 

b { 

float newangle = j graph- >GetEdge (choice, endjunc) - 


:3>Get Angle (choice) 
siewangle) ; 

5 


} 


float angle = ComputeAngle (* (endAngle [seedID] . end () -2) , 
deltaEnergy - = angleW * angle * angle; 

} 

changeAngle = true; 


seedID) ) 


//If the edge is no longer occupied after removal, increase 
// energy. 

if ( (edgeOccupation [j graph- >GetEdgeID (choice, endjunc) ] .size () ==1) && 
(*edgeOccupation [j graph- >GetEdgeID (choice, endjunc) ] .begin ( ) == 

deltaEnergy += unusedW * edgeLength; 


if ( (deltaEnergy < O.Of) || (exp (-deltaEnergy/ ( tempera tureConst* 
temperature) ) ) > GETRAND) 

{ 

// Remove the end junction. 
ASSERT (! paths [seedID] .empty ()) ; 

edgeOccupation [j graph- >GetEdgeID (choice, endjunc)] . erase (seedID) ; 
paths [seedID] .pop_back() ; 

if (changeAngle) 

{ 

ASSERT (endAngle. size () > 0) ; 
endAngle [seedID] .pop_back () ; 
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} 

// if the chosen junction was the hypo/rad separation for 
// the seedling, unmark the separation, 
if (RHseparation [seedID] == choice) 

{ 

RHseparation [seedID] = -1; 

} 

if (RHseparation2 [seedID] == choice) 

{ 

RHseparation2 [seedID] = -1; 

} 

pathlength [seedID] - = edgeLength; 
energy += deltaEnergy; 

} 

Q // Add edge if the selected junction is not the end junction itself 

=,p else if ((choice != endjunc) && (edgeOccupation [j graph- >GetEdgeID (endjunc, 

Lqhoice)] . find (seedID) == edgeOccupation [j graph- >GetEdgeID (endjunc , choice)] . end())) 

£ { 

// Add edge 


// * > 


* - 


; b j // start end choice 

;b ASSERT (j graph- >IsEdge (endjunc, choice)) ; 

,™ bool newseparation = false; // Whether choice adds a hypo/rad 

^Separation or not. 

:: :f bool separationcomplete = false; 

'™ float deltaEnergy; 

; * 3 float edgeLength = j graph- >GetEdge (endjunc, choice) ->m_length; 

// If the edge is not already occupied, the energy goes down, 
if (edgeOccupation [j graph- >GetEdgeID (endjunc, choice)] .empty () ) 
deltaEnergy - = unusedW * edgeLength; 

int numneighbors = 0; 

// See if choice is a separation point by checking for a 
// neighbor edge that is short (and has a degree one at the other 

for(int j=0; j < jgraph- >GetSize ( ) ; j++) 

{ 

if ( j graph- >IsEdge (choice, j ) ) 

{ 

numneighbors ++ ; 

} 

} 

// Mark if this is a new hypo/rad separation, 
if (RHseparation [seedID] == -1) 

{ 

if ( (numneighbors > 2) && (pathlength [seedID] > 20. Of)) 

{ 

newseparation = true; 


end??) 
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} 

} 

else if ( (RHseparation2 [seedID] == -1) && (pathlength [seedID] - 
separationlen [seedID] > 20. Of)) 

{ 

deltaEnergy -= separationW; 
separationcomplete = true; 

} 

if (edgeLength > minEdgeLength) 

{ 

if ( ! endAngle [seedID] . empty ( ) ) 

{ 

float newangle = j graph - >Ge t Edge (end junc, choice) - 

>GetAngle (endjunc) ; 

float angle = ComputeAngle (endAngle [seedID] . front () , 

newangle) ; 

deltaEnergy += angleW * angle * angle; 
£3 // egraph->GetEdge ( jgraph- 

:,|GetEdgeID (paths [seedID] [paths [seedID] .size()-2] , endjunc) , j graph ->GetEdge ID (endjunc , 
i Choice) ) ; 

tf» } 

M= changeAngle = true; 

,p } 

; s y if ( (deltaEnergy < O.Of) || (exp ( -deltaEnergy/ (temperatureConst* 

temperature) ) ) > GETRAND) 

a < 

m II mark that this seedling has occupied the edge 

edgeOccupation [j graph- >GetEdgeID (endjunc, 

Choice) ] . insert (seedID) ; 

jig ASSERT (edgeOccupation [j graph- >GetEdgeID (endjunc , choice) ] . size ( ) 

.1^= start junc . size ()) ; 

paths [seedID] .push_back (choice) ; 

if (changeAngle) 

endAngle [seedID] ,push_back (j graph- >GetEdge (endjunc, 

choice) ->GetAngle (choice) ) ; 

pathlength [seedID] += edgeLength; 
if (newseparation) 

{ 

RHseparation [seedID] = choice; 
separationlen [seedID] = pathlength [seedID] ; 

} 

if (separationcomplete) 

RHseparation2 [seedID] = choice; 

energy += deltaEnergy; 

} 

} 

// Adjust temperature, 
temperature *= 0.99f; 
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} // Simulated annealing loop, 
return paths; 

} 

float SeedlingAnalyzer : : ComputeAngle (float lastAngle, float newAngle) 

{ 

// Compute the angle between lastAngle and newAngle. 

float angle = fabs (newAngle - lastAngle); 
if (angle > PI) 

angle = DPI - angle; 

return (PI -angle) ; 

} 


us 

:. : 3 

rU 
□ 

■•n 

-j. « 

□ 
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// File: Image. cpp 

// Author: Yusaku Sako 

// Date: 07/11/99 

#include "stdafx.h" // for Windows 
#include " Image . h" 

#include "ImageThinning.h" // for thinning a binary image 
#include "Connectivity-Graph. h" 
#include "Const .h" 
#include <math.h> 

ttifndef ROUND 

#define ROUND (x) ( (int) (x+0 . 5) ) 
#endif 



Exlmage* Exlmage :: Crop (int lx, int ly, int rx, int ry) 

{ 

Image* dupimage = duplicate_Image (m_image) ; 
Q Image* newimage = :: crop (dupimage, ly, lx, ry-ly+1, rx-lx+1) ; 


f\ Exlmage* returnimage = new Exlmage; 

returnimage->m_image = newimage; 


..g return returnimage; 

a 

ru 

^Exlmage* Exlmage :: Trans formToHSV( int numlevels) 

^ Exlmage* newimg = (Exlmage*) Get Copy () ; 

£q float norm [3] ; 

; ia norm[0] = numlevels; normfl] = numlevels; norm [2] = numlevels; 

I'a newimg- >m_image = colorxform (newimg- >m_image, HSV, norm, NULL, 1) ; 


} 


return newimg; 


Exlmage* Exlmage : :Get2DHistogram (int xcolor, int ycolor, int size) 

{ 

Exlmage* histlmg = new Exlmage; 

histlmg->m_image = new_Image (PGM, GRAYSCALE, 1, size, 
Size, CVIP_BYTE, REAL) ; 

const int HIGH = 255; 
const int BIAS = 12 8; 

unsigned char ** xpix = (unsigned char**) m_image->image — ptr [xcolor] ->rptr; 

unsigned char ** ypix = (unsigned char**) m_image->image — ptr [ycolor] ->rptr; 

unsigned char ** hpix = (unsigned char**) histlmg ->m_image->image_ ptr [0] ->rptr; 

int kx, ky, hx, hy, max, kk; 

// Clear image. 

for(int i=0; i<size; i++) 

for(int j=0; j<size; j++) 

{ 

hpix[i] [j] = (unsigned char)0; 
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} 

max = 0; 
kx = 1; 
ky = l; 

for(i=0; i<size; i++) 

{ 

for(int j=0; j<size; j++) 
{ 

hy = (HIGH - ypixfi] [j])/ky; 
hx = (xpixfi] [j] ) /kx; 
if (hpix[hy] [hx] < HIGH) 

hpix [hy] [hx] ++; 
if (max < hpix[hy] [hx] ) 

max = hpix [hy] [hx] ; 

} 

} 

for(i=0; i<size; i++) 

{ 

i'3 for(int j=0; j<size; j++) 

£ { 

Lfj if (hpix [i] [j] !=0) 

4 { 

|,4 kk = hpixfi] [j] *HIGH/max+BIAS; 

[ a E if(kk > HIGH) 

\jj hpixfi] [j] = (unsigned char) HIGH; 

fy else 

hpix[i] [j] = (unsigned char)kk; 

} 

t "| return histlmg; 

:= ExImage* Exlmage :: Threshold (unsigned char thresh) 

{ 

Exlmage *returnlmage = new Exlmage; 

returnlmage->m_image = threshold_segment (m_image, thresh, CVIP_N0) ; 

/* 

<inputlmage> - pointer to Image structure 
<threshval> - threshold value 
<thresh_inbyte> 

- CVIP_N0 apply threshval directly to image data; 

- CVIP_YES threshval is CVIP_BYTE range; remap to image 

data range before thresholding. 

*/ 

return returnlmage; 

} 

Exlmage* Exlmage :: Threshold (int rlow, int rhigh, int glow, int ghigh, int blow, int 
bhigh) 

{ 

Exlmage* threshlmg = new Exlmage; 

threshlmg->m_image = new_Image (PGM, GRAYJSCALE, 1, getNoOf Rows_Image (m_image) , 
getNoOf Cols_Image (m_image) , CVIP_BYTE, REAL) ; 


n } 
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unsigned char** origR = (unsigned char**) m_image-> image _j?tr [0] ->rptr; 
unsigned char** origG = (unsigned char**) m_image-> image j>tr [1] ->rptr; 
unsigned char** origB = (unsigned char**) m_image->image_ptr [2] ->rptr; 

unsigned char** dest = (unsigned char**) threshlmg->m_image->image_ptr [0] ->rptr; 

for(int y=0; y<getNoOf Rows_Image (m_image) ; y++) 

for(int x=0; x<getNoOf Cols__ Image (m_image) ; x++) 

{ 

bool flag = true; 

if (origR [y] [x] < rlow | | origR [y] [x] > rhigh) 

flag = false; 
if (origG [y] [x] < glow | | origG [y] [x] > ghigh) 

flag = false; 
if (origB[y] [x] < blow | | origB [y] [x] > bhigh) 

flag = false; 

if (flag) 

n dest [y] [x] = 255; 

''5 else 

?i dest [y] [x] = 0; 


i 

return threshlmg; 


: ExImage* Exlmage : : Threshold (ThreshParams thresh) 

;f return Threshold (thresh . rmin, thresh. rmax, thresh. gmin, thresh. gmax, thresh. bmin, 

4?hresh .bmax) ; 

I 

^P/ Perform Thresholding segmentation based on histogram 
lllxlmage* Exlmage : :HistogramThreshold ( ) 

{ 

Exlmage *returnlmage = new Exlmage; 
returnlmage->m_image = hist_thresh_gray (m_image) ; 
return returnlmage; 

} 

Exlmage* Exlmage: :HistogramEqualization () 

{ 

Exlmage * returnlmage = new Exlmage; 
Image* temp; 

temp = histeq (m_image , 0) ; 
//temp = histeq(temp, 1) ; 
//temp = histeq(temp, 2) ; 
returnlmage ->m_image = temp; 

returnlmage ->m_image = remap__Image (temp, CVIP_BYTE, 0, 255); 


} 


return returnlmage; 


Exlmage* Exlmage : :EdgeDetect (int type) 
{ 
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Exlmage * re turn Image = new Exlmage; 

returnlmage->m_image = edge_detect_setup (m_image, 1) ; 
return returnlmage; 

} 

// Perform skeletonization on the image. 
Exlmage* Exlmage :: Thin () 

{ 

Exlmage * returnlmage = new Exlmage; 

returnlmage ->m_image = : : ImageThinning (m_image) ; 

return returnlmage; 

} 

Exlmage* Exlmage : :MorphOpen (int kerneltype, int kernelsize, 
int kernelheight # int kernelwidth) 

{ 

Exlmage * returnlmage = new Exlmage; 

returnlmage- >m_image = : :MorphOpen (m_image , kerneltype, kernelsize, kernelheight, 
kernelwidth) ; 
□ // DEBUG 

! I5 CString msg; 

; s fj msg. Format ("Image params: width=%d height=%d bands=%d colorspace=%d" , returnlmage- 
:.gm_image->image — ptr [0] ->cols, returnlmage->m_image->image_ptr [0] ->rows, returnlmage- 
i*^m_image - >bands , returnlmage - >m_image - >color_space ) ; 

//AfxMessageBox (msg) ; 
•/'i return returnlmage; 

ft 

^Exlmage* Exlmage :: Mo rphClose (int kerneltype, int kernelsize, 
! !~ int kernelheight, int kernelwidth) 

Exlmage * returnlmage = new Exlmage; 
;-Q returnlmage ->m_image = : :MorphClose (m_image, kerneltype, kernelsize, kernelheight, 

^kernelwidth) ; 

return returnlmage; 

} 

Exlmage* Exlmage : :MorphDi late (int kerneltype, int kernelsize, 
int kernelheight, int kernelwidth) 

{ 

Exlmage *returnlmage = new Exlmage; 

returnlmage ->m_image = : :MorphDilate (m_image, kerneltype, kernelsize, kernelheight, 
kernelwidth) ; 

return returnlmage; 

} 

Exlmage* Exlmage :: MorphErode (int kerneltype, int kernelsize, 
int kernelheight, int kernelwidth) 

{ 

Exlmage *returnlmage = new Exlmage; 

returnlmage ->m_image = :: MorphErode (m_i mage, kerneltype, kernelsize, kernelheight, 
kernelwidth) ; 

return returnlmage; 

} 

// Convert a color image to gray scale based on luminance. 
Exlmage* Exlmage :: ConvertToGray scale (int maxvalue) 
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{ 

Exlmage *returnlmage = new Exlmage; 

returnlmage ->m_image = :: CVIP luminance (m_image , maxvalue, 
CVIP_YES, CVIP_NO) ; 
return returnlmage; 

} 

// Perform median filter on the image. 
Exlmage* Exlmage : : MedianFilter (int kernelsize) 

{ 

Exlmage *returnlmage = new Exlmage; 

returnlmage ->m_image = : :median__f ilter (m_image, kernelsize); 
return returnlmage; 

} 

Exlmage* Exlmage :: CutOut (int lx, int ly, int rx, int ry) 
{ 

ASSERT (lx>=0 ) ; 

ASSERT (ly>=0) ; 
□ ASSERT (rx<getNoOf Cols_Image (m_image) ) ; 
i s Q ASSERT (ry<getNoOf Rows_Image (m_image) ) ; 

m 

!.p Image* dupimg = duplicate_Image (m_image) ; 

: ;e s, unsigned char** pix = (unsigned char**) dupimg- >image_ ptr [0] ->rptr; 

:y for (int y=ly; y<=ry; y++) 

ru { 

for (int x=lx; x<-rx; x++) 

{ 

pix[y] [x] = 0U; 

} 

} 


□ 


} 


Exlmage *returnlmage = new Exlmage; 
returnlmage ->m_image = dupimg; 
return returnlmage; 


Exlmage* Exlmage :: Minus (Exlmage* inimage, int lx, int ly, int rx, int ry) 

{ 

Image* diffimg = duplicate_Image (m_image) ; 

unsigned char** pixl = (unsigned char**) dif f img->image_ptr [0] ->rptr; 
unsigned char** pix2 = (unsigned char**) inimage- >m_image- >image_j?tr [0] ->rptr ; 

for (int y=ly; y<=ry; y++) 

{ 

for (int x=lx; x<=rx; x++) 

{ 

int diff = pixl [y] [x] - pix2 [y-ly] [x-lx] ; 
diff = (diff<0) ? 0 : diff; 
pixl [y] [x] = diff ; 

} 

} 

Exlmage* returnlmage = new Exlmage; 
returnlmage- >m_image = diffimg; 
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} 


return returnlmage; 


Exlmage* Exlmage :: Blend (Exlmage* inimage, float weight_orig, float weight_in, bool 
invert_orig, bool invert_in) 

{ 

// Image sizes must match. 

ASSERT (getNoOf Cols_Image (m_image) == getNoOf Col s_Image (inimage ->m_image) ) ; 
ASSERT (getNoOf Rows_Image (m_image) == getNoOf Rows_Image (inimage ->m_image) ) ; 
ASSERT ( (weight_orig >^ 0.0) && ( weight_orig <= 1.0)); 

Image* newimage = new__Image (PPM, RGB, 3, getNoOf Rows_Image (m_image) , 
getNoOf Col s_Image (m_image) , CVIP_BYTE, REAL) ; 


int numbands = getNoOf Bands_Image (inimage ->m_image) ; 

char **srclR = m_image- >image_ptr [0] - >rptr ; 

char **src2R = inimage - >m_image - >image_ptr [0] ->rptr; 

unsigned char **destR = (unsigned char**) newimage- >image_ptr [0] ->rptr; 
. c3 q char **srclG = m_image-> image jptr [1] ->rptr ; 

char **src2G = (numbands>l) ? inimage - >m_image - >image_ptr [1] ->rptr 
fsm_image->image — ptr [0] ->rptr; 

s = g unsigned char **destG = (unsigned char**) newimage- >image_ptr [1] ->rptr; 

/"T char **srclB = m_image->image_ptr [2] ->rptr ,* 

char **src2B = (numbands>2) ? inimage ->m_image->image_ptr [2] ->rptr 
~>m_image->image_ptr [0] ->rptr; 

^ unsigned char **destB = (unsigned char**) newimage- >image_ ptr [2] ->rptr; 

< U 

unsigned char srclr, src2r, srclg, src2g, srclb, src2b; 

for (int y=0; y<getNoOf Rows_Image (m_image) ; y++) 

for (int x=0; x<getNoOf Cols_Image (m_image) ; x++) 

{ 

if (invert_orig) 


a 

: : 

3 


srclr = -srclR[y] [x] -1; 
srclg = -srclG[y] [x]-l; 
srclb = -srclB[y] [x]-l; 


{ 


} 

else 

{ 

srclr = srclR[y] [x] ; 
srclg = srclG[y] [x] ; 
srclb = srclB [y] [x] ; 

} 

if (invert_in) 

{ 

src2r = -src2R[y] [x]-l; 
src2g = -src2G[y] [x] -1; 
src2b = -src2B[y] [x] -1; 

} 

else 

{ 

src2r = src2R [y] [x] ; 

src2g = src2G [y] [x] ; 

src2b = src2B [y] [x] ; 


inimage - 


inimage - 
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destRty] [x] 

char) ( (weight_orig* (float) (srclr) +weight_in* (float) src2r) /2) 
destG[y] [x] 

char) ( (weight_orig* (float) (srclg) +weight_in* (float) src2g) /2 ) 
destB [y] [x] 

char) ( (weight_orig* (float) (srclb) +weight_in* (float) src2b) /2) 

} 


(unsigned 
(unsigned 
(unsigned 


Exlmage *ret = new Exlmage; 
ret->m_image = newimage; 
return ret; 


Exlmage* Exlmage :: Mask (Exlmage* maskimage, ColorType bgcolor) 

{ 

Exlmage * outimage = (Exlmage*) Get Copy () ; 


unsigned char** mpix = maskimage- >GetPixelArray (0) ; 

unsigned char** srcR = GetPixelArray (0) 

q unsigned char** srcG = GetPixelArray (1) 

i.p unsigned char** srcB = GetPixelArray (2 ) 

Lf| unsigned char** dstR = outimage- >GetPixelArray (0) 

; ,p unsigned char** dstG = outimage- >GetPixelArray (1) 

|.I unsigned char** dstB = outimage- >GetPixelArray (2) 


1*3 

i. u 


HO 
Q 

□ 


for(int y=0; y<maskimage->GetHeight () ; y++) 
{ 

for(int x=0; x<maskimage- >GetWidth ( ) ; x++) 


{ 


if (mpix [y] [x] > 0) 

{ 


} 

else 
{ 


dstRfy] [x] 
dstG[y] [x] 
dstB[y] [x] 


dstR[y] [x] 
dstG[y] [x] 
dstB [y] [x] 


srcR [y] [x] 
srcG[y] [x] 
srcB [y] [x] 


bgcolor . r; 
bgcolor .g; 
bgcolor .b; 


} 

return outimage; 


ObjectList label_Objects2 (Image *imageP, Image **labelP, unsigned background); 


// Return objects found in the image. 
BinaryObjects* Exlmage: :GetObjects () 

{ 

Image* label Image; 
ObjectList objlist; 

objlist = label_Objects2 (m_image, &labellmage, 0) ; 
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// Make sure objlist is not NULL 
if ( ! objlist) 

objlist = new_LL(); 
BinaryObjects *binobjs = new BinaryObjects (objlist, label Image ) ; 
return binobjs; 

} 

int Exlmage : : FindPeakValue (int band) 

{ 

// Find peak. 

unsigned char **pix = (unsigned char**)m_image->image — ptr [band] ->rptr; 
long* histogram = new long [256] ; 

for (int i=0; i<256; i++) 
histogram [i] = 0; 

for(int y=0; y<GetHeight ( ) ; y++) 

{ 

for(int x=0; x<GetWidth ( ) ; x++) 
histogram [pix [y] [x]]++; 

m } 


// Find the peak. 
\\j int high = 0; 

?fi for(i=0; i<256; i++) 

;. y ( 

r=i if (histogram [i] > histogram [high] ) 

™ high = i; 


delete [] histogram; 

□ 

^ return (high) ; 

} 

ColorType Exlmage : : GetAverageColor ( ) 

{ 

long rsum = 0, gsum = 0, bsum = 0; 

unsigned char** pixR = (unsigned char**) m_image->image — ptr [0] ->rptr ; 
unsigned char** pixG = (unsigned char**) m_image->image_ptr [1] ->rptr ; 
unsigned char** pixB = (unsigned char**) m_image-> image _ptr [2] ->rptr; 

for (int y=0; y<getNoOf Rows_Image (m_image) ; y++) 

for (int x=0; x<getNoOf Cols_Image (m_image) ; x++) 

{ 

rsum += pixR [y] [x] ; 
gsum += pixG [y] [x] ; 
bsum += pixB [y] [x] ; 

} 

ColorType color; 

color. r = rsum/ (getNoOf Rows — Image (m_image) *getNoOf Cols_Image (m_image) ) ; 
color. g = gsum/ (getNoOf Rows_Image (m_image) *getNoOf Cols_Image (m_image) ) ; 
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} 


color .b = bsum/ (getNoOf Rows_Image (m_image) *getNoOf Cols_Image (m_image) ) ; 
return color; 


// Dump properties of all objects, 
void BinaryOb j ec t s : : Dump ( ) 

{ 

head_LL (m_objectlist) ; // set linked list pointer to the head 
for(int i=0; i<size_LL (m_objectlist) ; i++) 

{ 

Object *obj = ((Object *) retrieve_LL (m_obj ectlist) ) ; 
CString msg; 

msg. Format ( "label=%d; R=%d G=%d B=%d; xmin=%d ymin=%d xmax=%d ymax=%d;\n 
eigenratio=%f ; orientation=%f ; horizontal cog=%f vertical cog=%f; area=%f" , 

obj->label, obj ->pixel .r, obj ->pixel .g, obj ->pixel . b, obj->x_min, ob j - 
>y_min, obj->x_max, obj->y_jnax, obj ->prop . eig_ ratio, obj ->prop .orientation, ob j - 
>prop . h_cog , obj - >prop . v_cog , obj - >prop . area ) ; 

AfxMessageBox (msg) ; 
l ^ next_LL (m_obj ectlist) ; 

hD } 

II 

: =Int _cdecl Matchlnt (void* content, void* lookforP) 
'"'J return (*( (int*) content) == *( (int*) lookforP) ) ; 

Remove objects whose sizes are equal to or less than minarea. 
r^oid BinaryObjects : :FilterOut (int minarea) 

& 

if (m_objectlist->listlength == 0) 
return; 

getProp_Objects (m_objectlist , m_labellmage) ; 

head_LL (m_obj ectlist) ; 
previous_LL (m_obj ectlist) ; 

//int** pix = (int**)m_labellmage->imagej>tr [0] ->rptr; 


for ( ; ; ) 
{ 


>label) 


Object* obj = ((Object *) retrieveNext_LL (m_obj ectlist) ) ; 
if (obj ->prop . area < minarea) 

{ 

// Erase pixels belonging to this object from the label image 
for (int y=obj ->y__min; y<=obj - >y_max; y++) 

for (int x=obj ->x_min; x<=obj ->x_max; x++) 

{ 

if ( ( (int*) m_labellmage->image — ptr [0] ->rptr [y] ) [x] == obj 

{ 

( (int*)m_labellmage->image_ptr [0] ->rptr [y] ) [x] = 0; 
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} 

removenext_LL (m_objectlist ) ; 

} 

else 

next_LL (m_obj ectlist) ; 

if (istail_LL (m_obj ectlist) ) 
break; 

} 

} 

Exlmage* BinaryOb j ect s : : CreateContour Image ( ) 

{ 

Image* returnimg; 

returnimg = new_Image ( PBM , BINARY, 1, getNoOf Rows__Image (m_label Image ) , 
getNoOfCols_Image (m_label Image ) , CVIPJBYTE, REAL) ; 

/* 

: ,~ returnimg - new_Image (TIF, BINARY, 

; s - 1, getNoOfRows_Image (m_label Image) , getNoOf Cols_Image (m_label Image) , 
!^ CVIP BYTE, REAL) ; 

m */ ~ 

;"T // Extract chain code from each object, and draw them onto returnimg. 
// dump labelled image 


CSt ring msg; 
/* 


£3 msg. Format ("# bands=%d width=%d height=%d format=%d type=%d, space=%d'\ 

: : fti_label Image ->bands , 

□ getNoOf Col s_Image (m_labellmage) , getNoOf Rows_Image (m_label Image) , 

•jn_l abel Image ->image_format , 

q m_labellmage- > image _jptr [0] ->data_jtype, m_label Image ->color_ space) ; 

r'3 AfxMessageBox (msg) ; 

msg. Format ("# bands=%d width=%d height=%d format=%d type=%d, space=%d", returnimg- 
>bands , 

getNoOf Cols_Image (returnimg) , getNoOf Rows_Image (returnimg) , returnimg- 

>image_f ormat , 

returnimg- >image_ptr [0] ->data_type, returnimg- >color_space) ; 
AfxMessageBox (msg) ; 
*/ 

head_LL (m_obj ectlist) ; 

for( ;; ) 
{ 

int ray_x; 

Object *obj = ( (Object*) retrieve_LL (m_obj ectlist) ) ; 

//shootRay (m_label Image, obj->label, &ray_x, &ray__y, obj->x_min, obj->y__min, 
ob j - >x_max , ob j - >y_max) ; 

for (int c=obj ->x_min; c<=obj ->x_max; C++) 

{ 

if ( ( (int*) m_labellmage->imagejptr [0] ->rptr [obj ->y_min] ) [c] == ob j - 

>label) 

{ 
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ray_x = c ; 
break; 

} 

} 

ChainCode* cc = new_ChainCode (obj ->y_min / ray_x, obj->label); 

if (build_LineChainCode (cc, m_ labellmage, obj->x_min, obj->y_min, obj ->x_max, ob j - 
>y_max) ==0) 

AfxMessageBox ( "Error building chain code!"); 

draw_ChainCode (cc, returnimg) ; 
if (istail_LL (m_obj ectlist) ) 

break; 
next_LL (m_obj ectlist) ; 

} 

Exlmage* newimg = new Exlmage; 
newimg->m_image = returnimg; 

return newimg; 

} 

;ExImage* BinaryOb j ects : : Createlmage ( ) 
„ Image* returnimg; 

^ returnimg = new_Image (PGM, GRAY_SCALE, 1, getNoOf Rows_Image (m_label Image ) , 

getNoOfCols_Image (m_label Image) , CVIP_BYTE, REAL) ; 


for(int y=0; y<getNoOf Rows_Image (m_label Image ) ; y++) 

for(int x=0; x<getNoOf Cols_Image (m_label Image) ; x++) 

■ { 

if ( ( (int*) m_labellmage->image_ptr [0] ->rptr [y] ) [x] > 0) 
y] returnimg- >image_ptr [0] ->rptr [y] [x] = 255U; 

□ else 

uQ returnimg- >image_ptr [0] ->rptr ty] [x] = 0; 

Q } 


} 


Exlmage* newimg = new Exlmage; 
newimg ->m_image = returnimg; 

return newimg; 


Connect ivityGraph<Edge *>* BinaryOb j ects : : Create JunctionGraph (Object Info* obj info, 
std : : vector<Junction*>& junclist /* output */) 

{ 

// Find all junctions. 

junclist = ComputeJunctions (obj info) ; 

TRACE ( "There are %d pre- junctions\n" , junclist . size ()) ; 

ConnectivityGraph<Edge *>* eg = new Connect ivityGraph<Edge *> (false, 
junclist . size ( ) ) ; 

//int *degreelist = new int [junclist . size ()] ; 

// Copy degrees. 

//for(int i=0; junclist . size () ; i++) 

//{ 
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// 

//} 


degreelist [i] = junclist [i] ->degree; 


/ / Connect nodes . 

for(int i = 0; i< junclist . size () ; i++) 

{ 

int c = 0 ; 
float oldlen; 

while (junclist [i] ->m_neighborPixels . size () != 0) 

{ 

Edge *edge = new Edge; 

TRACE ("j [%d] 's next 
>m_neighborPixels .begin () ) ) ; 

int nextjunc = 

>m_neighborPixels .begin ()) , i, *edge) ; 

// If nextjunc is i itself, then this is a terminal loop. 

// Simply ignore it. 

if (nextjunc !- i) 


neighbor is %d\n", i 
FindNext June t ion ( j unc list, 


* (junclist [i] 

* (junclist [i] 


'-'=3 

:.§dge . 


{ 


ru 


en- 


TRACE ("Inserting edge (%d, %d)\n", i, nextjunc); 
// If there is already an edge (i, nextjunc) 
if ( !cg->IsEdge (i, nextjunc) ) 

{ 

oldlen = edge->m_length; 
cg->InsertEdge (i, nextjunc, edge) ; 

} 

else if (edge->m_length < oldlen) 


ignore the new 


{ 


} 


cg->DeleteEdge (i, nextjunc) ; 
cg->InsertEdge (i, nextjunc, edge) ; 


>GetNeighbor (i) , 
>GetNeighbor (nextjunc) ) 

} 

else 

{ 

neighbor 


// Erase the neighbor that leads to nextjunc and 

// erase the neighbor that leads to i so as to avoid 

// duplicate processing of the same edge. 

TRACE ("Erasing neighbor %d from j [%d] and %d from j[%d]\n" , edge- 
edge- >GetNeighbor (next j unc) , nextjunc) ; 

junclist [i] ->m_neighborPixels. erase (edge->GetNeighbor (i) ) ; 
junclist [nextjunc] ->m_neighbor Pixels . erase (edge- 


//degreelist [i] ; 

C + +; 


// Erase the neighbor pixel that leads to the loop 

// In case of loop, edge->GetNeighbor (0) returns the first 

// and edge->GetNeighbor (1) returns the last neighbor 
junclist [i] ->m_neighborPixels .erase (edge->GetNeighbor (0) ) ; 
junclist [i] ->m_neighborPixels . erase (edge- >GetNeighbor (1) ) ; 
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} 


TRACE ("There are %d post edges\n", cg->GetNumEdges () ) ; 
return eg,* 


int BinaryObjects : :FindNext Junction (std: :vector< Junction *> junclist, int neighbor, int 

start, Edge& edge) 

{ 

// FOR DEBUG 

FILE* out = fopen(" junction.txt", "w") ; 

int** pix = (int**) m_labellmage->imagej>tr [0] ->rptr; 

int deg; 

int x = junclist [start] ->m_x; 
int y = junclist [start] ->m_y; 

edge. m_j unci = start; 

// Keep track of all coordinates, 
edge . m_xarray . clear ( ) ; 
^pj edge .m_yarray. clear () ; 

edge . m_xarray . push__back (x) ; 
;: P edge.m yarray.push back (y) ; 

fU edge .m_xmin = 1000000; 

n edge . m_xmax = - 1 ; 

q edge.m_ymin = 1000000; 

;j| edge.m_ymax = -1; 


II FOR DEBUG 

fprintf(out, "finding the next junction for junction %d\n", start); 


do 
{ 


fprintf (out, "neighbor=%d\n" , neighbor) ; 
switch (neighbor) 

{ 

case 0 : 

x = x+1; 

y = y; 

edge .m_length += l.Of; 
break; 
case 1: 

x = x+1; 
y = y+i; 

edge.m_length += 1.41421356f; 
break; 
case 2 : 

x = x; 

y = y+i; 

edge . m_length += 1 . 0 f ; 
break; 
case 3 : 

X = X-l; 


22727/04060 


13 


Listing 2 


y = y+i; 

edge . m_length += 1.41421356f; 
break; 
case 4 : 

x - x-1 ; 

y = y; 

edge .m_length += l.Of; 
break; 
case 5 : 

x = x-1; 

y = y-i; 

edge.m_length += 1.41421356f; 
break; 
case 6 : 

x = x; 

y = y-i; 

edge . m_length += l.Of; 
break; 
case 7 : 

x = x+1; 

:.n y = y-i; 

T~ edge.m_length += 1.41421356f; 

* a L; break; 

;^ } 

; 3 if (x < edge . m_xmin) 

edge.m_xmin = x; 
; - if (x > edge . m_xmax) 

edge . m_xmax = x; 
l P if (y < edge ,m_ymin) 

J™ edge . m_ymin - y; 

: -=^ if (y > edge .m_ymax) 

Vy edge . m_ymax = y ; 

O edge . m_xarray . push_back (x) ; 

edge ,m_yarray .pushjback (y) ; 

//Is the current position a junction? 

// neighbor 
// 5 6 7 
// 4 0 
// 3 2 1 

int nextNeighbor; 

deg=0; 

if ( (pix[y] [x+1] != 0) && (neighbor != 4)) 
{ 

deg++; 

nextNeighbor = 0; 

} 

if ( (pix[y+l] [x] != 0) && (neighbor != 6)) 
{ 

deg++ ; 

nextNeighbor = 2 ; 
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in 


f ( (pix[y] [x-1] != 0) && (neighbor != 0)) 
deg++; 

nextNeighbor = 4; 
f ( (pix[y-l] [x] != 0) && (neighbor != 2)) 
deg++; 

nextNeighbor = 6; 

f ( (pix[y+l] [x+1] != 0) &&(pix[y] [x+1] ==0) &&(pix[y+l] [x] ==0) && (neighbor != 5)) 
deg++; 

nextNeighbor = 1; 

f ( (pix[y+l] [x-1] != 0) &&(pix[y+l] [x] ==0) && (pix [y] [x-1] ==0)&& (neighbor != 7)) 
deg++; 

nextNeighbor = 3; 

f ( (pix[y-l] [x-1] != 0) ScSc (pix [y-1] [x]==0)&& (pix [y] [x-1] ==0)&& (neighbor != 1) ) 
deg++; 

nextNeighbor - 5; 

f ( (pix [y-1] [x+1] != 0) && (pix [y-1] [x] ==0) && (pix [y] [x+1] ==0) && (neighbor != 3)) 
deg++; 

nextNeighbor = 7; 


fprintf(out, "degree=%d nextneighbor=%d\n" , deg, nextNeighbor); 
neighbor = nextNeighbor; 
} while (deg == 1) ; 

// Compute the angle the edge forms at each junction (with respect to x-axis) 

// 1. angle at junction 1 

int chainlen = edge . m_xarray . size () ; 

int dx = edge .m_xarray [min (4 , chainlen-1) ] - edge .m_xarray [0] ; 
int dy = edge .m_yarray [min (4 , chainlen-1)] - edge ,m_yar ray [0] ; 

if (dx == 0) 
{ 

if (dy < 0) 

edge ,m_anglel = HPI; 

else 

edge ,m_anglel = -HPI; 


} 

else 


edge .m_anglel = atan2 ( (double) -dy, (double) dx) ; 


dx = edge. m_xarray [max (chainlen- 5, 0)] - edge .m_xarray [chainlen-1] ; 
dy = edge ,m_yarray [max (chainlen-5 , 0)] - edge .m_yarray [chainlen-1] ; 


if (dx == 0) 
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{ 

if (dy < 0) 

edge ,m_angle2 = HPI; 

else 

edge .m_angle2 = -HPI; 

} 

else 

edge .m_angle2 = atan2 ( (double) -dy, (double)dx); 

// Return the ID of the junction 
for(int i=0; i< junclist . size ( ) ; i++) 

{ 

if ( (junclist [i] ->m_x = = x) && (junclist [i] ->m_y == y) ) 

{ 

edge . m_j unc2 = i ; 
f close (out) ; 
return i; 

} 

} 

,rj for(i=0; i<junclist . size ( ) ; i++) 

5 { 

fprintf(out, "junc %d=%d %d\n" / i, junclist [i] ->m_x, junclist [i] ->m_y) ; 

} 

f close (out) ; 
ASSERT (FALSE) ; 
return -1; // Error! ! ! 


';SonnectivityGraph<f loat>* BinaryObjects : : CreateEdgeGraph (ConnectivityGraph<Edge *>* jg, 
!:-QbjectInf o* objinfo) 

// Edge Graph 

□ Connect ivityGraph<f loat>* eg = new ConnectivityGraph<float> (false, jg- 

CiGetNumEdges () ) ; 

//CString msg; 

//msg. Format ("junction graph: number of edges=%d", jg->GetNumEdges ( ) ) ; 
//AfxMessageBox (msg) ; 

// For each edge-to-edge connection, compute the angle between the edges, 
int JGsize = jg->GetSize () ; 
for(int i=0; i<JGsize; i++) 

{ 

for (int j=i+l; j< JGsize; j++) 

{ 

if (jg->IsEdge (i, j)) 

{ 

for (int k=0; k< JGsize; k++) 

{ 

if (jg->IsEdge(j ( k) && (i ! = k) ) 
{ 

// Compute the angle between i,j,k 
// i 
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k))) 


jg->GetEdge ( j , k) ->GetAngle ( j ) ; 


// 

// / \ 
// / \ 
// / \ 
// j k 

if ( ! eg->IsEdge ( jg->GetEdgeID (i, j), jg->GetEdgeID ( j , 


{ 


>GetEdgeID ( j , k) , angle) ; 


float angle = jg->GetEdge (i , j ) ->GetAngle ( j ) 

if (angle < O.Of) 

angle += DPI; 
if (angle > PI) 

angle = DPI - angle; 
eg->InsertEdge (jg->GetEdgeID (i, j) , jg- 

//CString msg; 

//msg. Format ("Inserting edge %d %d for %d %d 


jg->GetEdgeID(i, j) , jg->GetEdgeID ( j , k) , i, j , k) ; 

//Af xMessageBox (msg) ; 


•JO 


Q 
Q 


k)) 


jg->GetEdge (i,k) ->GetAngle (i) 


>GetEdgeID (i, k) , angle); 


} 
} 

else if (jg->IsEdge (i, k) && (j != k) ) 
{ 

// Compute the angle between i,j,k 

// j 
// 

// / \ 

// / \ 
// / \ 
// i k 

if ( ! eg->IsEdge ( jg->GetEdgeID (i, j), jg->GetEdgeID (i , 
{ 

float angle = jg->GetEdge (i , j ) ->GetAngle (i) 

if (angle < O.Of) 

angle += DPI; 
if (angle > PI) 

angle = DPI - angle; 
eg->InsertEdge ( jg->GetEdgeID (i, j) , jg- 


CString msg; 

msg. Format ("Inserting edge %d %d for %d %d %d M , 
jg->GetEdgeID(i, j) , jg->GetEdgeID (i , k) , i, j, k) ; 

Af xMessageBox (msg) ; 

} 

} 

} 


} 

return eg; 
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} 

Exlmage* BinaryObjects : : Create June t ionlmage (std : : vector<Obj ectlnf o*> objinfo) 

{ 

Image* returnimg; 

returnimg = new_Image (PPM, RGB, 3, getNoOf Rows_Image (m_label Image ) , 
getNoOfCols_Image (m_label Image) , CVIP_BYTE, REAL) ; 

head_LL (m_objectlist) ; 

int** srepix = (int**) m_labellmage->image j>tr [0] ->rptr; 

unsigned char** destpixR = (unsigned char**) returnimg- > image j>tr [0] ->rptr; 
unsigned char** destpixG = (unsigned char**) returnimg- >image_ptr [1] ->rptr; 
unsigned char** destpixB = (unsigned char**) returnimg- >image_ptr [2] ->rptr; 

for(int i=0; i<size_LL (m_obj ectlist) / i++) 

{ 

Object* obj = (Object*) retrieve_LL (m_obj ectlist ) ; 
for (int y = obj->y_min; y <- obj->y_max; y++) 

{ 

for (int x = obj->x_min; x <= obj->x_max; x++) 

{ 

if (srepix [y] [x] == obj->label) 

{ 

destpixR [y] [x] = 255; 
destpixG [y] [x] = 255; 
destpixB [y] [x] = 255; 

} 

} 

} 

linked_list* jlist = obj info [i] -> junctions ; 
head_LL (jlist) ; 
; r0 for(int j=0; j <size_LL ( j list) ; j++) 

S < 

: :=^ Junction* junc = ( Junction*) retrieve_LL (j list) ; 

for (int k=0; k< junc- >m_degree; k++) 

{ 

int jx = junc->m_x, jy = junc->m_y; 
switch ( junc->m_degree) 

{ 

case 3 : 


o 

■■J3 


a 


case 4 


case 5 


destpixR [ jy] 

[jx] 


255; 

destpixG [ jy] 

[jx] 


0; 

destpixB [jy] 

[jx] 


0; 

break; 




destpixR [jy] 

[jx] 


255; 

destpixG [ jy] 

[jx] 


255; 

destpixB [ jy] 

[jx] 


0; 

break; 




destpixR [ jy] 

[jx] 


0; 

destpixG [ jy] 

[jx] 


255; 

destpixB [ jy] 

[jx] 


0; 

break; 




destpixR [ jy] 

[jx] 


0; 
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destpixGtjy] [jx] = 255; 
destpixB [ jy] [ jx] = 255; 
break; 
case 7 : 

destpixR [ jy] [jx] = 0; 
destpixG [jy] [ jx] = 0; 
destpixB [jy] [jx] = 255; 
break; 

} 

} 

next_LL ( j list) ; 

} 

next_LL (m_obj ectlist) ; 

} 

Exlmage *newimg = new Exlmage; 
newimg->m_image = returnimg; 

return newimg; 

i 

in 

^Hoid Ob j ect Info : : Dump ( ) 
:: r CStrmg msg; 

I.U msg. Format ( "xmin=%d xmax=%d ymin=%d ymax=%d\nperimeter=%d area=%d eigenratio=%f 

fBrientation=%f \nxcenter=%f ycenter=%f numjunc=%d" # 

;s xmin, xmax, ymin, ymax, perimeter, area, eigenratio, orientation, xcenter, 

Reenter, numjunc) ; 

AfxMessageBox (msg) ; 

l=^td: : vector<Obj ect Info* > BinaryObjects : :GetObjectInf o (bool buildChaincode) 

std: : vector<0bjectlnfo*> infolist; 

if (m_objectlist->listlength == 0) 
return infolist; 

getProp_Objects (m_obj ectlist , m_labellmage) ; 

head_LL (m_objectlist) ; 

for(int i=0; i<size_ LL (m_obj ectlist) ; i++) 

{ 

int ray_x ; 

Object* obj = ( (Object*) retrieve_LL (m_obj ectlist ) ) ; 

//shootRay (m_label Image, obj->label, &ray_x, &ray_y, obj->x_min, obj->y_min, 
obj - >x_max , obj - >y_ max) ; 

for (int c=obj ->x_min; c<=obj - >x_max; C++) 

{ 

if ( ( (int*)m_labellmage->image — ptr [0] ->rptr [obj ->y_min] ) [c] == ob j - 

>label) 

{ 

ray_x = c ; 
break; 
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} 

/* 

for(int y=obj ->y_min; y<=obj ->y_max; y++) 
{ 

for(int x=obj ->x_min; x<=obj ->x_max; x++) 

{ 

CString msg; 

msg. Format ( "i=%d objlabel=%d, label=%d xmin=%d xmax=%d ymin=%d 
ymax=%d xstart=%d ystart=%d" / i, obj->label, ( (int*) m_l ab el Image ->image_j?tr [0] - 
>rptr [y] ) [x] , obj ->x_min, obj ->x_max, obj ->y_min, obj ->y_max) ; 

AfxMessageBox (msg) ; 

} 

} 

*/ 

Objectlnfo* info = new Objectlnfo; 
:«=s if (buildChaincode) 

3 < 

i'^ info->chain = new_ChainCode (obj ->y_min, ray_x, obj->label); 

^ if ( !build_LineChainCode (inf o- >chain, m_label Image, obj ->x_min, obj - 

; : ?y_min , obj - >x_max , obj - >y_max ) ) 

2 AfxMessageBox ( "Error building chain code"); 

*' r , //print_ChainCode (inf olist [i] ->chain, "chaininfo.txt"); 

,y inf o->perimeter = inf o->chain->no_of_vectors ; 

} 

info->area = obj ->prop . area; 
U? info->label = obj->label; 

P inf o->eigenratio = obj ->prop . eig_ratio; 

inf o->orientation = obj ->prop .orientation ; 
; jP inf o->xcenter = obj ->prop . h__cog; 

O inf o->ycenter = obj ->prop . v_cog; 

q info->xmin = obj->x_min; 

inf o- >xmax = obj - >x_max ; 
info->ymin = obj->y_min; 
inf o- >ymax = obj - >y_max ; 

//inf o->junctions = Compute June t ions (inf o) ; 
inf olist .push_back(info) ; 
next_LL (m_objectlist) ; 

} 

return inf olist; 

std: : vector<ExImage*> BinaryObj ects : : CreateColorlmages ( std : : vector<0bjectlnf o*> infolist, 
ImagelO* original) 


ru 


std: : vector<ExImage*> imglist; 

int** lab = (int* *)m_label Image ->image_j>tr [0] ->rptr; 
unsigned char** origR = original - >GetPixelArray (0) 
unsigned char** origG = original - >GetPixelArray (1) 
unsigned char** origB = original ->GetPixelArray (2) 

for(int i=0; i<inf olist . size () ; i++) 
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Objectlnfo* info = infolist [i] ; 

// Create a color image for the object. 

Image* returnimg = new_Image (PPM, RGB, 3, info->ymax - info->ymin + 1, 
info->xmax - info->xmin + 1, CVIP_BYTE, REAL) ; 

unsigned char** destR = (unsigned char**) returnimg- >image_ptr [0] ->rptr; 
unsigned char** destG = (unsigned char**) returnimg- >image_ptr [1] ->rptr; 
unsigned char** destB = (unsigned char**) returnimg- >image_j?tr [2] ->rptr; 


int count = 0; 


for(int y=info->ymin; y<=inf o- >ymax; y++) 

{ 

for (int x=inf o->xmin; x<=info->xmax; x++) 

{ 

int ny = y - info->ymin; 
int nx = x - info->xmin; 


□ 


□ 

i!3 


} 

} 

Exlmage* 
newimage - 


if (labty] [x] == info->label) 
{ 

destR [ny] [nx] = origR [y] [x] 
destG [ny] [nx] = origG [y] [x] 
destB [ny] [nx] = origB [y] [x] 
count++; 


} 

else 

{ 


newimage 
>m_image 


destR [ny] [nx] 
destG [ny] [nx] 
destB [ny] [nx] 


= new Exlmage; 
= returnimg; 


imglist .push_back (newimage) ; 


} 

return imglist; 


(unsigned char)0; 
(unsigned char)0; 
(unsigned char)0; 


Exlmage* BinaryObjects : : CreateColorlmage (Objectlnfo* info, ImagelO* original) 

{ 

int** lab = (int**) m_l abe 1 Image -> image jptr [0] ->rptr; 
unsigned char** origR = original- >GetPixelArray (0) ; 
unsigned char** origG = original - >GetPixelArray (1) ; 
unsigned char** origB = original ->GetPixelArray (2 ) ; 

// Create a color image for the object. 

Image* returnimg = new_Image (PPM, RGB, 3, info->ymax - info->ymin + 1, 
info->xmax - info->xmin + 1, CVIP_BYTE, REAL) ; 

unsigned char** destR = (unsigned char**) returnimg- >image_j?tr [0] ->rptr; 
unsigned char** destG = (unsigned char**) returnimg- >image_ptr [1] ->rptr; 
unsigned char** destB = (unsigned char**) returnimg- >image_ptr [2] ->rptr; 
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int count = 0; 

for(int y=inf o->ymin; y<=inf o->ymax; y++) 

{ 

for (int x=inf o->xmin; x<=info->xmax; x++) 

{ 

int ny = y - info->ymin; 
int nx = x - info->xmin; 

if(lab[y][x] info->label) 
{ 

destR[ny] [nx] = origR [y] [x] ; 
destG [ny] [nx] - origG [y] [x] ; 
destB [ny] [nx] = origB [y] [x] ; 
count++; 

} 

else 

{ 

destR[ny] [nx] = (unsigned char)0; 

q destG [ny] [nx] = (unsigned char)0; 

. 5 ~ destB [ny] [nx] = (unsigned char)0; 


} 


} 

Exlmage* newimage = new Exlmage; 
newimage ->m_image = returnimg; 

return newimage; 


I'ii 
1 

Mbftd: : vector<Junction*> BinaryObj ects : : Compute Junctions (Obj ectlnfo* objinfo) 

i 

std: :vector<Junction*> junclist; 


: Pa 


Q int** pix = (int**) m_labellmage->image_ptr [0] ->rptr ; 

int label = obj info- >label ; 

for (int y=obj inf o->ymin; y<=obj inf o- >ymax; y++) 

for (int x=obj inf o->xmin; x<=objinfo->xmax; x++) 

{ 

std: : set<int> neighbor; 

if (pix[y] [x] == label) 
{ 

// neighbor 
// 5 6 7 
// 4 0 
// 3 2 1 

int deg=0; 

if ( (x!=objinfo->xmax) && (pix[y] [x+1] == label)) 
{ 

deg++; 

neighbor . insert ( 0 ) ; 

} 

if ( (y ! =objinfo->ymin) && (pix[y-l][x] == label)) 
{ 
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deg++; 

neighbor . insert (6) ; 

} 

if ( (x!=objinfo->xmin) && (pix[y][x-l] == label)) 

{ 

deg++; 

neighbor . insert (4 ) ; 

} 

if ( (y!=objinfo->ymax) && (pix[y+l][x] == label)) 

{ 

deg++; 

neighbor . insert (2) ; 

} 

if { (y I =obj info- >ymin) && (x!=objinfo->xmax) && (pix [y-1] [x+1] 

label) ) 

{ 

if ( (neighbor. find (6) ==neighbor . end ( ) ) 
(neighbor . find (0) ==neighbor . end ( ) ) ) 

{ 

deg++; 

l'r* neighbor . insert ( 7 ) ; 

s , » 

/- if ( (y ! =obj info->ymin) && (x!=objinfo->xmin) && (pix [y-1] [x-1] 

! tabel) ) 

if ( (neighbor, find (4) ==neighbor . end ( ) ) 
^neighbor. find (6) ==neighbor . end ( ) ) ) 

% { 

'J deg++; 

uO neighbor . insert ( 5 ) ; 

□ } 
.0 } 

i"3 if ( (y ! =obj info->ymax) && (x ! =obj inf o- >xmin) (pix [y+1] [x-1] 

rlabel) ) 

{ 

if ( (neighbor . find (2 ) —neighbor . end ( ) ) 
(neighbor . find (4 ) ==neighbor . end ( ) ) ) 

{ 

deg++; 

neighbor . insert (3) ; 

} 

} 

if ( (y i =obj info->ymax) && (x ! =obj inf o->xmax) && (pix [y+1] [x+1] 
label)) . 

{ 

if ( (neighbor . find (0) ==neighbor . end ( ) ) 
(neighbor . find (2) —neighbor . end () ) ) 

{ 

deg++; 

neighbor . insert ( 1 ) ; 

} 

} 

// Modify neighborhood list based on the following rules: 
// if the neighbor is odd and there is an adjacent neighbor, 
// remove the neighbor from the list. 
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/* 

std :: set<int> :: iterator it = neighbor . begin () ; 
int c = 0 ; 

while (it != neighbor . end () ) 

{ 

if(*it % 2 == 1) 
{ 

int n = *it; 

int nl = (n != 0) ? n-1 : 7; 
int n2 = (n != 7) ? n+1 : 0; 
if ( (neighbor. find (nl) != 
(neighbor . find (n2 ) != neighbor . end ()) ) 

{ 

it = neighbor . erase (it) ; 
deg--; 

} 

} 

else 

it + +; 

^ if(c>10) 

AfxMessageBox ( "ERROR" ) ; 

C++; 


neighbor . end ( ) ) 


ii 


m 


ru 

Q 

a 


} 

return j unci is t; 


if (deg != 2) 
{ 


Junction *junc = new Junction; 
junc->m_degree = deg; 
j unc - >m_x = x ; 
junc->m_y = y; 

j unc ->m_neighbor Pixels = neighbor; 
junclist .push_back ( junc) ; 


void BinaryObjects : :Merge Junctions (ConnectivityGraph<Edge *> *jgraph, float 
me rgeDi stance) 

{ 

// juncs keeps track of junctions that belong to "short" edges. 

std: : set<int> juncs; 

std : : set<int > : : iterator it ; 

Point* juncpos = new Point [ jgraph->GetSize ()] ; 
for(int i=0; i< jgraph->GetSize ( ) ; i++) 

{ 

// Find edges that are short enough and put their junctions into juncs. 
// Remove those edges. 

for(int j=i+l; j < jgraph->GetSize ( ) ; j++) 

{ 

if (j graph- >IsEdge (i , j ) ) 

{ 
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# 


if ( j graph- >GetEdge (i # j ) ->m_length < mergeDistance) 

{ 

juncpos[i] = j graph ->Get Edge (i, j ) ->GetEndPosition (i) ; 

juncs . insert (i) ; 

juncs . insert ( j ) ; 

j graph- >Delet eEdge ( i , j ) ; 

} 


a 

S3 


std: :vector<std: :vector<int> > groups; 

// Look at junctions' connectivity and divide them into groups 
//of connected junctions, 
while ( ! juncs . empty ( ) ) 

{ 

std: :vector<int> list; 
groups .push_back (list) ; 

int current = * juncs . begin () ; 
groups. rbeginO ->push_back (current) ; 
juncs. erase (juncs. begin () ) ; 


^ std: :vector<int> worklist; 

=P do 

S { 

iy it = juncs .begin () ; 

while (it != juncs. end () ) 

y { 

Cn if ( jgraph->IsEdge (current, *it) ) 

□ { 

;,Q worklist .push_back(*it) ; 

it = juncs . erase (it) ; 

else 

it++; 

} 

if (worklist . empty ( ) ) 
break; 

current = *worklist . begin () ; 
worklist . erase (worklist . begin () ) ; 
groups . rbegin () ->push_back (current) ; 
} while ( ! worklist . empty ( ) ) ; 

} 

int* newjunclndex = new int [groups . size ()] ; 
Point* newjuncPos = new Point [groups . size ()] ; 

// Compute the centroid for each group and make it a new junction. 
for(i=0; i<groups . size ( ) ; i++) 

{ 

float sumx = O.Of; 
float sumy = O.Of; 

for (int j=0; j <groups [i] . size ( ) ; j++) 

{ 

sumx += juncpos [groups [i] [j]].x; 
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sumy += juncpos [groups [i] [j] ] .y; 

} 

sumx /= groups [i] . size () ; 

sumy /= groups [i] . size () ; 

Point pt (ROUND (sumx) , ROUND (sumy) ) ; 

newjuncPos [i] = pt; 

newj unc Index [i] = groups [i] [0]; 

} 

for(i=0; i<groups . size () ; i++) 

{ 

for(int j=0; j <groups [i] . size ( ) ; j++) 

{ 

// Modify all edges that are connected to groups [i] [ j ] 
for(int k=0; k< jgraph- >GetSize () ; k++) 

{ 

if ( jgraph- >IsEdge (groups [i] [ j ] , k) ) 

{ 

Edge* edge = j graph- >GetEdge (groups [i] [j ] , k) ; 
//AddSegmentEnd (edge, newjuncPos [i] , newjunclndex [i] , k, 


QergeDi stance) ; 

m 

Aspect to x-axis) 


// Compute the angle the edge forms at each junction (with 

// 1. angle at junction 1 

int chainlen = edge ->m_xar ray . size ( ) ; 

int endlndex; 


if (newjunclndex [i] == edge ->m_j unci) 
^ endlndex = 0; 

else 

endlndex = chainlen - 1; 


int dx = newjuncPos [i] .x - edge ->m_xarray [endlndex] ; 
int dy = newjuncPos [i] .y - edge ->m_yar ray [endlndex] / 

if (dx == 0) 
{ 

if(dy > 0) 
{ 

if (newjunclndex [i] •== edge ->m_j unci) 
edge->m_anglel = HPI; 

else 

edge->m_angle2 = HPI; 


} 

else 

{ 


} 

else 

{ 


if (newjunclndex [i] == edge ->m_j unci) 
edge->m_anglel = -HPI; 

else 

edge->m_angle2 = -HPI; 


i f ( newj unc Index [ i ] = = edge - >m_ j unc 1 ) 


22727/04060 


26 


Listing 2 


(double) dx) ; 


(double) dx) ; 

} 

} 

} 

} 

} 

void BinaryObjects : : Adds egment End (Edge* edge, Point newPos, int newlndex, int otherlndex, 
float mergeDistance) 

{ 

// Connect the edge to the newly formed junction. 

ASSERT ( edge ->m_j unci == otherlndex || edge->m_junc2 == otherlndex); 

if (edge- >m_j unci == otherlndex) 

r==: edge->m junc2 = newlndex; 

else 

^ edge ->m_j unci = newlndex; 

= n 

iM/ Helper function for GetObj ectlnf o . 

avoid BinaryObjects: : Compute Junctions (Ob j ectlnf o* objinfo, Junction **junctions, int 
L&jnum j unc ) 

! ? 

Q int xmin = obj info->xmin; 

int ymin = obj info- >ymin; 
Q int xmax = obj info- >xmax; 

i;3 int ymax = obj info- >ymax; 

// Convert the chain code into a list of (x,y) . 
int *xP=NULL, *yP=NULL; 

getXY_ChainCode (obj info- >chain, &xP, &yP) ; 

// count is a counter for each coordinate in the chain code, 
int** count = new int* [ymax-ymin+1] ; 
for (int i=0; i<ymax-ymin+l; i++) 

{ 

count [i] = new int [xmax-xmin+1] ; 
// Reset counter, 
for (int j=0; j <xmax-xmin+l ; j++) 
count [i] [j ] =0; 

} 

numjunc = 0; 

for (i = 0 ; i<obj inf o->chain->no_of_vectors; i++) 

{ 

// Increment numjunc if the junction degree is greater than 2. 
if (count [yP [i] -ymin] [xP [i] -xmin] == 2) 

{ 

numjunc++; 


edge->m_anglel = atan2 ( (double) -dy, 

else 

edge->m_angle2 = atan2 ( (double) -dy, 
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} 

(count [yP [i] -ymin] [xP [i] -xmin] ) ++; 


// No junction. 
if(numjunc == 0) 

{ 

* junctions = NULL; 
return; 

} 


Junction* newjunctions = new Junction [numjunc] ; 


int jc = 0; 

for(i=0; i<ymax-ymin+l; i++) 

for (int j=0; j <xmax-xmin+l ; j++) 

{ 

if (count [i] [j] >= 3) 
{ 

} «j newjunctions [jc] .degree = count [i] [j]; 

newjunctions [jc] .x = j+xmin; 
newjunctions [jc] .y = i+ymin; 
jc++; 

} 

} 

* junctions = newjunctions; 
ASSERT ( j c == numjunc); 


in 


ru 


II Clean up. 
delete [] xP; 
0 delete [] yP; 


;Bypedef HashTable *ObjectHash; 


#define HASH SIZE 251U 


#def ine hash_Object (label) ( ( (unsigned) (label) ) %HASH_SIZE) 

static void addto_Objects (ObjectHash hashP, int next^abel, Color pixel, int y_ pos, int 
x_jpos) ; 

static void update_Obj ects (ObjectHash hashP, int object_label , int y_j?os, int x_pos) ; 
static void combine_Objects (Obj ectHash hashP, int b, int c, int *xmin, int *xmax, int 
*ymin, int *ymax) ; 

static ObjectList hash2List_0bj ects (ObjectHash hashP) ; 
static void makegraymap_Obj ects (ColorHistogram *chP) ; 

static Matrix *color2index_Image (ColorHistogram *chP, Image *imageP ); 
/* 

* label recycling routines 
*/ 

static void initLabelStack (void) ; 
static int getNextLabel (void) ; 
static void recycleLabel (int label); 
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/* 

* global variables used by label recycling routines 
*/ 

static int label_count; 
static Stack label_stackP; 

ObjectList 
label_0bjects2 ( 

Image *imageP # 
Image **labelP, 
unsigned background 
) 

{ 

register int x, y=0, i, j; 
int A, B, C, D, rows, cols; 

byte *pixarray2, *pixarrayl; /* pixel arrays */ 

int *labarrayl, *labarray2, *rowP; /* label arrays */ 

int xmin, xmax, ymin, ymax, next_label; 

unsigned *A_LABEL, B_LABEL, C_LABEL, D_LABEL, MAX_LABEL; 
q ObjectList listP; 

ObjectHash hashP; 
fa Object *objP; 

ColorHistogram *chP; 
{.2 ColorHistObject *mapP; 

r Matrix *matrixP; 

ffj ROI *roiP; 
j!y const char *fn = "label"; 

chP = new_ColorHist () ; 

;L. S if (getNoOfBands_Image (imageP) > 1) { 

u ^ compute_ColorHist (chP, imageP, 256); 

^ matrixP = color2index_Image (chP, imageP); 


{ 

makegraymap_Objects (chP) ; 
matrixP = getBand_JEmage (imageP, 0) ; 


mapP = chP->histogram; 
if (mapP==NULL) return NULL; 

rows = getNoOfRows_Image (imageP) ; 
cols = getNoOf Cols_Image (imageP) ; 

*labelP = new_Image (PGM, GRAY_SCALE, 1, rows, cols, CVIP_INTEGER, REAL); 

pixarray2 = (unsigned char*) getRow_Matrix (matrixP, 0) ; 
labarray2 = (int *) getRow_Image (* label P, 0, 0) ; 

initLabelStack () ; 

hashP = new_HT(HASH_SIZE) ; 

/* 

* handle special case of first row 
*/ 



else 
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for(x=0; x < cols; x++) 

if( pixarray2 [x] != background ) { 

if( (x==0) || (pixarray2 [x-1] != pixarray2 [x] ) ) { 
next_label = getNextLabel () ; 

addto_Objects (hashP, next_label, mapP [pixarray2 [x] ] .pixel, y, x) ; 
labarray2 [x] = next_label; 

} 

else { 

update_Objects (hashP, labarray2 [x-1] , y, x) ; 
labarray2 [x] = labarray2 [x-1] ; 

} 

} 

for(y=0; y < rows-1; y++) 

{ 

pixarrayl = pixarray2 ,- 

pixarray2 = (unsigned char* )getRow_Matrix (mat rixP, y+1) ; 


Q labarrayl = labarray2 ; 

v3 labarray2 = (int *) getRow_Image (*labelP / y+1, 0) ; 

m 

ip for(x=0; x < cols-1; x++) 

M { 

A = pixarray2 [x+1] ; 

B = pixarray2 [x] ; 

fU C = pixarrayl [x+1] ; 


: - 1 


C3 


sis* 

n 


D = pixarrayl [x] ; 

A_LABEL = (unsigned *) &labarray2 [x+1] ; 
B_LABEL = labarray2 [x] ; 
C_LABEL = labarrayl [x+1] / 
D_LABEL = labarrayl [x] ; 

if ( X == 0) 
{ 

if (B != background) 
if (D == B) 
{ 


} 

else 

{ 


mapP[B] .pixel, y+1, x) 


} 

} // x==0 

if (A != background) 
{ 

if(D != A) 
{ 

if (B != A) 
{ 


B_LABEL = labarray2 [x] - D_LABEL; 
update_Objects (hashP, B_LABEL, y+1, x) ; 


next_label = getNextLabel () / 

addto_Objects (hashP, next_label, 
B_LABEL = labarray2 [x] = next_label; 
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next_label , mapP [A] .pixel, y+1, x+1) 


if (C != A) 
{ 

next_label = getNextLabel () ; 
addto_Objects (hashP, 

*A_LABEL = next_label; 
} // C!=A 
else 

{ 


X+l) ; 


*A_LABEL = C_LABEL; 
update_Objects (hashP, 

} // C==A 
} // B!=A 
else 

{ 

if (C == A) 
{ 

if (B LABEL == C LABEL) 


*A LABEL, 


r.3 

in 


:: !C_LABEL, &xmin, &xmax, 

ru 
a 
a 

j-^min, xmax-xmin+1, 

□ 

i++) { 

getRow_ROI (roiP, i, 0) ; 
getNoOfCols_ROI (roiP) ; j++, rowP++) 


{ 
} 

else 
{ 


*A_LABEL = B_LABEL / 

combine_Objects (hashP, 
&ymin, &ymax) ; 


B_LABEL, 


MAX_LABEL = MAX (B_LABEL, C_LABEL) ; 
* A_LABEL = MIN ( B_LABEL , C_LABEL ) ; 

roiP = new_ROI(); 

asgnImage_ROI (roiP, *labelP, xmin, 
ymax-ymin+1) ; 

for(i = 0; i < getNoOf Rows_ROI (roiP) / 
rowP = (int *) 

for(j=0; j < 


} 


if (*rowP == MAX_LABEL) 
*rowP = *A LABEL; 


delete_ROI (roiP) ; 
} // B_LABEL ! =C_LABEL 
} // C==A 
else 

*A_LABEL = B_LABEL; 

update_Objects (hashP, *A_LABEL, y+1, x+l); 
} // B!=A 
} // D!=A 
else 

{ 
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*A_LABEL = D_LABEL; 

update_Objects (hashP, *A_LABEL, y+1, x+1) ; 
} // D==A 
} // A! =background 

else // addition by Yusaku Sako. 7/26/99 

{ 

if (B != background) 

if((B == C) && (D != B) ScSc (B_LABEL != C_LABEL) ) 
{ 

combine_Objects (hashP, B_LABEL , C_LABEL, &xmin, 

&xmax, 

&ymin, &ymax) ; 


xmin+1 , 

□ 

:.H| 

M?owP++) 

=.5 


MAX_LABEL = MAX (B_LABEL , C_JLABEL) ; 

B_LABEL = labarray2 [x] = MIN (B_LABEL, C__LABEL) ; 

roiP = new_ROI(); 

asgnImage_ROI (roiP, *labelP, xmin, ymin, xmax- 

ymax-ymin+1) ; 

for(i=0; i < getNoOf Rows_ROI (roiP) ; i++) { 

rowP = (int *) getRow_ROI (roiP, i, 0) ; 
for(j=0; j < getNoOf Cols_ROI (roiP) ; 


} 


if (*rowP == MAX_LABEL) 
*rowP = B LABEL; 


: =erase 
■>3 


delete_ROI (roiP) ; 

update_Objects (hashP, B_LABEL, y+1, 


x) 


} 


} // A— background 
} // for x 
} // for y 


// 


if (hashP->table [hashP->key] ) 

listP = hash2List_0bjects (hashP) ; 

else 

listP = NULL; 


/* 

fprintf (stderr, "Displaying Object Statistics ... \n\n" ) ; 

fprintf (stderr, "Number of objects found = %u.\n\n" / size_LL (listP) ) ; 

print_Objects (listP) ; 

fprintf (stderr, "\n\n Press the < ENTER > key to continue"); 

getchar () ; 

*/ 


// Bug fix by YS (3/24/00) . 

// Don't wanna delete matrixP when number of bands is 1 
if (getNoOf Bands_Image (imageP) > 1) { 
delete^ Matrix (matrixP) ; 

} 


return listP; 
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} 

static ObjectList hash2List__Objects (ObjectHash hashP) 

{ 

register int i, first = 0; 
ObjectList listP; 

fort;;) { 

if (hashP->table [first] ) break; 
f irst++; 

} 

listP = hashP->table [first] ; 

f or (i=f irst+1; i < size_HT (hashP) ; i++) 
if (hashP->table [i] ) { 

* (listP->tailP) = * (hashP->table [i] ->headP->nextP) ; 

listP->tailP = hashP->table [i] ->tailP; 

listP->listlength += size__LL (hashP->table [i] ) ; 
□ delete_Link (hashP->table [i] ->headP->nextP) ; 

;>S delete_Link (hashP->table [i] ->headP) ; 

:.fj /* delete_Link (hashP->table [i] ->headP->nextP) ; */ 

v3 free (hashP->table [i] ) ; 


i 

return listP; 


^tatic void addto_Objects (ObjectHash hashP, int next_label, Color pixel, int y_j?os, int 

setKey_HT (hashP, hash_Object (next_label) ) ; 
^3 addObject_HT ( hashP, new_Object (next_label, pixel, y_j?os, x_ pos) ); 

static void update_Objects ( ObjectHash hashP, int object_label, int y_pos, int x_j?os) 

{ 

Object *obj ; 

const char *fn = "update^ Objects" ; 

setKey_HT(hashP,hash_Object (object_label) ) ; 

if (f indObject_HT (hashP, match_Ob j ect , &object_label) ) { 
obj = (Object*)getObject_HT(hashP) ; 

obj -> x_min = MIN(obj -> x_min, xjjos) ; 

obj -> x_max = MAX (obj -> x_max, xjpos) ; 

obj -> y_min = MIN(obj -> y_min, y_j>os) ; 

obj -> y_max = MAX (obj -> y_max, y_j?os) ; 

} 


static void combine_Objects (ObjectHash hashP, int b, int c, int *xmin, int *xmax, int 
*ymin, int *ymax) 
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{ 

dlink *linkP; 
Object *bP, *cP; 

if (b < c) 
{ 

int temp; 
temp = b; 
b = c; 
c = temp; 

} 

setKey_HT(hashP,hash_Object (b) ) ; 
f indNextObject_HT(hashP, match_Object , &b) ; 
bP = (Object* )getNextObject_HT(hashP) ; 
removeNextObj ect__HT (hashP) ; 

setKey_HT(hashP,hash_Object (c) ) ; 
f indObject_HT (hashP, match_Obj ect , &c) ; 
cP = (Object*)getObject_HT(hashP) ; 

□ 

*xmin = bP->x_min; 
|f| *xmax = bP->x_max; 

*ymin = bP->y_min; 
M *ymax = bP->y_max; 

\j\ cP -> x_min = MIN(cP 

i=y cP -> x_max = MAX(cP 

' cP -> y_min = MIN(cP 

Q cP -> y_max = MAX(cP 

delete_Object (bP) ; 

£ recycleLabel (b) ; 

i j 

static void makegraymap__Obj ects ( ColorHistogram *chP ) 

{ 

register int i; 

chP->histogram = (ColorHistObject *) malloc (256*sizeof (ColorHistObj ect) ) ; 

for(i=0; i < 256; i++) 

assign_Color(chP->histogram[i] .pixel, i, i, i) ; 

} 


static Matrix *color2index_Image (ColorHistogram *chP, Image *imageP ) 

{ 

Matrix *matrixP; 

Color pixel; 

ColorHashTable *chtP; 

unsigned rows, cols; 

register int i; 

byte *mP, *rP, *gP, *bP; 


-> x_min, bP->x_min) ; 

-> x_max, bP->x_max) ; 

-> y_min, bP->y_min) ; 

-> y_max, bP->y__max) ; 
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if (chP->no_of_colors > 256) { 
return NULL; 

} 

rows = getNoOf Rows_Image (imageP) ; 
cols = getNoOf Cols_Image (imageP) ; 

rP = (unsigned char*) get Row_Image (imageP, 0, RED); 
gP = (unsigned char*) getRow_Image (imageP, 0, GRN) ; 
bP = (unsigned char*) getRow_lmage (imageP, 0, BLU) ; 

matrixP = new_Matrix (rows , cols, CVIP_BYTE, REAL); 
mP = (unsigned char*) getRow_Matrix (matrixP, 0) ; 

chtP = hist2Hash_ColorHT(chP) ; 

for(i=0; i < rows*cols; i++, mP++, rP++, gP++, bP++) { 
assign_Color (pixel, *rP, *gP, *bP) ; 

if( (*mP = lookUpColor_ColorHT ( chtP, pixel )) == -1 ) 
j "? return NULL; 

:| } 

j =n delete_ColorHT(chtP) ; 
h= return matrixP; 

4 

Id 

ru 

-static void initLabelStack (void) 

rfj label_stackP = new_Stack(); 

: ==~ label_count=l ; 

3 

^tatic int getNextLabel (void) 

int next_label, *labelP; 

if (isempty_Stack(label_stackP) ) 

next_label = label_count++; 

else { 

labelP = (int*)pop_Stack(label_stackP) ; 
next_label = *labelP; 
free (labelP) ; 

} 

return next_label; 

} 

static void recycleLabel (int label) 

{ 

int *labelP = (int *) malloc (sizeof (int) ) ; 
*labelP = label; 

push_Stack (label_stackP, labelP) ; 

} 
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* 


// Connectedness. Used in ImageTh inning function, 
int cconc(int inb[9]) 

{ 

int icn = 0; 

for(int i=0; i<8; i+=2) 
{ 

if (inb[i]==0) 

if(inb[i+l] == 255 || inb[i+2] == 255) 
icn++; 

} 

return icn; 

} 

// Performs thinning of a binary image 
Image* ImageThinning (Image *img) 

{ 

// Copy original image 
Image *newimg = duplicate_Image (img) ; 
._, if(lnewimg) 

AfxMessageBoxC'Not enough memory!"); 
unsigned char **pix = (unsigned char **) newimg- > image j>tr [0] ->rptr; 
int m = 100; int ir = 1, ia [9] , ic [9] ; 


! PI 


m 


int numrows = img->image_ptr [0] ->rows; 
int numcols = img -> image _ptr [0] ->cols ; 


while (ir!=0) 

a { 

m ir=0; 


for (int iy=0; iy < numrows ; iy++) 

{ 

for (int ix=0; ix < numcols; ix++) 
{ 

if (pix[iy] [ix] != 255) 
continue; 

if (ix == numcols-1) 
ia[0] = 0; 

else 

ia[0] = pix[iy ] [ix+1] ; 
if (iy ==0 || ix == numcols-1) 
ia[l] = 0; 

else 

ia[l] = pix[iy-l] [ix+1] ; 
if (iy == 0) 

ia[2] = 0; 

else 

ia[2] = pix[iy-l] [ix ] ; 
if (iy ==o || ix == 0) 
ia[3] = 0; 

else 

ia[3] = pix[iy-l] [ix-1] ; 
if (ix == 0) 

ia[4] = 0; 
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else 

ia[4] = pix[iy ] [ix-1] ; 
if (iy == numrows - 1 | | ix == 0) 
ia[5] = 0; 

else 

ia[5] = pix[iy+l] [ix-1]; 
if (iy == numrows-1) 
ia[6] = 0; 

else 

ia[6] = pix[iy+l] [ix ]; 
if (iy == numrows-1 | | ix == numcols-1) 
ia[7] = 0; 

else 

ia[7] = pix[iy+l] [ix+1] ; 
for(int i=0; i < 8; i++) 

{ 

if (ia [i] == m) 

{ 

ia [i] = 255; ic [i] = 0; 

□ } 

hQ else 

in { 

:.p if(ia[i] < 255) 

ia[i] = 0; 

^ ic [i] = ia [i] ; 

} 

} // for 
ia [8] = ia [0] ; 
ic [8] = ic [0] ; 

if (ia [0] +ia [2] +ia [4] +ia [6] ==255*4) 
continue; 


: si 

ru 


int iv, iw; 

for(i=0, iv=0, iw=0; i<8; i++) 
{ 

if (ia[i]==255) 

iv++; 
if (ic [i] ==255) 

iw++; 

} 

if (iv<=l) 

continue; 
if (iw==0) 

continue; 
if (cconc (ia) ! =1) 

continue; 

int temppix - (iy ==0) ? 0 : pix[iy-l] [ix] ; 
if (temppix == m) 

{ 

ia[2] = 0; 

if (cconc (ia) ! = 1) 

continue; 
ia[2] = 255; 

} 

temppix = (ix =-0) ? 0 : pix[iy] [ix-1] ; 
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if (temppix == m) 


ia[4]=0; 

if (cconc (ia) ! =1) 


continue; 
ia[4] = 255; 


pix [iy] [ix] 
ir++; 
// for ix 


m; 


} // for iy 
m++; 
} // while 

for(int iy=0; iy<img->image_ptr [0] ->rows; iy++) 

{ 

for(int ix=0; ix<img->image_ptr [0] - >cols ; ix++) 


if (pix[iy] [ix] < 255) 
pixfiy] [ix] = 0; 


return newimg; 
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